A High-level Open Source Model for the Australian National Electricity Market (NEM)

Author

Grant Chalmers

Published

October 14, 2025

Introduction

The National Electricity Market (NEM), which spans roughly 40,000 km of transmission lines and cables, serves more than 23 million Australians across five interconnected regions - Queensland, New South Wales, Victoria, South Australia, and Tasmania.¹ As one of the world’s longest and most complex power systems, the NEM coordinates wholesale electricity generation and transports it from generators to large industrial users, local distributors, and ultimately to homes and businesses. This system relies on a real-time ‘pool’ market, with electricity output aggregated and scheduled every five minutes to efficiently balance demand and incentivise investment in fast-response technologies, such as batteries and gas peaking units.

As the NEM undergoes rapid transformation, driven by government policy promoting greater renewable energy integration and grid modernisation, advanced modelling tools have become essential for analysing the current state of the market and evaluating future scenarios.

This Quarto document presents the development of a high-level baseline model for the NEM power system using PyPSA (Python for Power System Analysis), an open-source framework for simulating and optimising modern energy systems. The model, referred to as High-level_NEM, is a 5-node model (one bus per region) based on actual 2024 generation and demand profiles at 30-minute resolution. A comparison of the 2024 baseline model with a future scenario is presented in this document. The future scenario has increased renewable energy and storage and associated reduction in thermal generation.

The workflow implemented in this Quarto document encompasses the following key steps:
- Data Import and Preprocessing: Loading and preparing network topology, generator, storage, time series data from actual NEM system records and relevant scenario input data.
- Network Construction: Building the PyPSA network by adding buses, loads, generators, storage units, and transmission lines.
- Temporal Resolution and Snapshots: Setting up the model to operate at a chosen temporal granularity, allowing for the study of daily and seasonal variability in demand and renewable generation.
- Optimisation and Scenario Analysis: Solving the baseline model using linear optimisation to determine the optimal dispatch of resources, and extending the analysis to a plausible future scenario with an aggressive increase in variable renewable energy and associated reduction in thermal generation.
- Results Visualisation and Verification: Generating plots and statistics to interpret model outputs, including network topology, dispatch profiles, generation mix, and importantly, curtailment.
- Scenario Functions: Providing reusable functions for projecting various scaled-up generation scenarios.

Research Objectives

Objective:

Develop a high-level 5-node baseline PyPSA model for the NEM and VRE/storage scale-up and thermal generation scale-down scenarios with:
1. One bus per NEM region with transmission interconnectors
2. 30-minute time resolution to capture some level of RE intermittency
3. Comparison of 2024 baseline model (0_2024_baseline) to a high penetration renewable energy scenario (8.3_VreCurtailReview)
4. Develop an interactive web prototype Streamlit app for visualising dispatch profiles.

Key Assumptions

Baseline 2024 model:
- 2024 demand (load)
- 2024 generation profiles
- 2024 renewable energy availability (capacity) factors (incl. hydro).
- Generator registered capacities sourced from Open Electricity facilities Inclusive of operating and committed, therefore includes Snowy 2.0, Hunter Gas and a significant increase in battery storage over what is currently operating.
- Time series data is at 30-minute resolution.
- Marginal price of wind, utility-scale solar, rooftop solar, battery storage and hydro² set to zero.
- Rooftop solar sourced from AEMO via nemosis python package (TYPE = ‘MEASUREMENT’).
- Efficiency of storage on the way in and out set at 0.917 (round trip 0.84) via 2024 ISP Inputs and Assumptions workbook.
- Starter battery (600MW Cellars Hill) added in TAS1 in order to be able to scale up (not currently operating, but has some level of state govt. support).
- Project EnergyConnect, Marinus stage 1 & VNI West implemented on top of existing interconnectors (constrained by nominal capacity).
- One generator of each major source and one battery setup per region/bus.

Note

Apart from the baseline model, 19 additional scenarios were created by varying the scale-up of VRE and storage, and the scale-down of thermal generation, ultimately culminating in the final 8.3_VreCurtailReview scenario.

By following this workflow, the document provides a high-level example for those learning power system analysis with PyPSA.

Important links/references:

TU Berlin: Data Science for Energy System Modelling
PyPSA Documentation and Components
PyPSA Earth Documentation
GitHub PyPSA Sources
PyPSA-PH: High-Resolution Open Source Power System Model for the Philippines
2024 Integrated System Plan (ISP)
Open Electricity
High‑Level NEM PyPSA Model Github Repository

01 Import

This section imports all necessary Python packages and libraries required for data processing, power system modelling, geospatial analysis, and visualisation. These packages provide the foundational tools for building, analysing, and visualising the High-level NEM baseline model and high penetration renewable energy model.

Show the Python Code

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
from matplotlib.dates import AutoDateLocator, ConciseDateFormatter
import pypsa
from pypsa.plot.maps.static import add_legend_patches
from datetime import timedelta
from datetime import datetime  
from dateutil.parser import parse
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from time import perf_counter
from great_tables import GT, md, system_fonts, nanoplot_options
from IPython.display import display, HTML

02 Create Network

This section covers the creation of the High-level NEM network, including loading network data, initialising the PyPSA network object, and preparing the structure for adding buses, loads, generators, storage units, and transmission lines (links & lines).

02.01 Load Network Data

Show the Python Code

# Load components
generators = pd.read_csv("data/generators.csv").set_index('name')
buses = pd.read_csv("data/buses.csv", index_col=0)
loads = pd.read_csv("data/loads.csv", index_col=0)
lines = pd.read_csv("data/lines.csv", index_col=0)
links = pd.read_csv("data/links.csv", index_col=0)
storage_units = pd.read_csv("data/storage_units.csv", index_col=0)

Show the Python Code

# Load time series data, which has 17,568 rows due to 2024 being a leap year (30-min resolution).
load_ts = pd.read_csv("data/loads_t.p_set.csv", index_col=0)
load_ts.index = pd.to_datetime(load_ts.index, errors="raise")
print(f'Length: {len(load_ts)}, Index data type: {load_ts.index.dtype}') 
# Below should be empty if all dates are valid
load_ts[~load_ts.index.notnull()]
generators_ts = pd.read_csv("data/generators_t.p_max_pu.csv", index_col=0)
generators_ts.index = pd.to_datetime(generators_ts.index, errors="raise")
print(f'Length: {len(generators_ts)}, Index data type: {generators_ts.index.dtype}') 


if not load_ts.index.equals(generators_ts.index):    
    raise ValueError("Time series indices are not aligned")

Length: 17568, Index data type: datetime64[ns]
Length: 17568, Index data type: datetime64[ns]

02.02 Initialise the network

Show the Python Code

n = pypsa.Network()

03 Add Network Components

This section details the process of adding key components to the network, including buses, loads, generators, storage units, and transmission lines. Each component is integrated with relevant time series data at 30-minute resolution.

Show the Python Code

# Add components
for name, row in buses.iterrows():
    n.add("Bus", name, **row.to_dict())

for name, row in loads.iterrows():
    n.add("Load", name, **row.to_dict())

for name, row in lines.iterrows():
    n.add("Line", name, **row.to_dict())

for name, row in generators.iterrows():
    n.add("Generator", name, **row.to_dict())

for name, row in links.iterrows():
    n.add("Link", name, **row.to_dict())
    
for name, row in storage_units.iterrows():
    n.add("StorageUnit", name, **row.to_dict())    


# Explicitly set index with correct frequency in a robust way (handles different temporal resolutions)
load_ts.index = pd.date_range('2024-01-01', periods=len(load_ts), freq='30min')
generators_ts.index = load_ts.index

# Explicitly set snapshots in PyPSA
n.set_snapshots(load_ts.index)

# Re-initialise snapshot_weightings explicitly
n.snapshot_weightings.loc[:, :] = n.snapshot_weightings.index.to_series().diff().dt.total_seconds().div(3600).fillna(0.5).values[:, None]

# Assign load and generator profiles again
n.loads_t.p_set = load_ts
n.generators_t.p_max_pu = generators_ts

print("Final snapshot weighting preview:")
print(n.snapshot_weightings.head())

print("Unique snapshot weights:", n.snapshot_weightings.iloc[:, 0].unique())

assert all(n.generators_t.p_max_pu.index == n.snapshots)

# Remove negative values in p_max_pu (Hydro is the culprit)
n.generators_t.p_max_pu = n.generators_t.p_max_pu.clip(lower=0.0, upper=1.0)

Final snapshot weighting preview:
                     objective  stores  generators
snapshot                                          
2024-01-01 00:00:00        0.5     0.5         0.5
2024-01-01 00:30:00        0.5     0.5         0.5
2024-01-01 01:00:00        0.5     0.5         0.5
2024-01-01 01:30:00        0.5     0.5         0.5
2024-01-01 02:00:00        0.5     0.5         0.5
Unique snapshot weights: [0.5]

Introduce a dummy ‘load shedding’ generator to help identify infeasible scenarios within the model. This approach is more informative than simply receiving a “not feasible” message from the solver, as it highlights which regions are affected. With this information, we can make targeted adjustments - such as increasing generation or storage in future iterations.

Set the marginal cost very high so it is used as a last resort.

Show the Python Code

# Add one unserved energy generator per load bus
# Acts like a dummy "load shedding" generator
for bus in n.loads.bus.unique():
    gen_name = f"{bus}-UNSERVED"
    n.add("Generator",
          name=gen_name,
          bus=bus,
          carrier="Unserved Energy",
          p_nom_extendable=True,
          p_nom=0,
          marginal_cost=10000,  # Very expensive fallback
          capital_cost=0,       # Optional: make it purely operational cost
    )

Show the Python Code

# Diagnostic output
print(f"Loaded {len(n.buses)} buses")
print(f"Loaded {len(n.loads)} loads with time series of shape {n.loads_t.p_set.shape}")
print(f"Loaded {len(n.generators)} generators with time series of shape {n.generators_t.p_max_pu.shape}")
print(f"Loaded {len(n.lines)} lines")
print(f"Loaded {len(n.links)} links")
print(f"Loaded {len(n.generators)} generators")
print(f"Loaded {len(n.storage_units)} storage units")

# Check total demand makes sense
print(f"Total demand:  {n.loads_t.p_set.sum(axis=1).max()}")

Loaded 5 buses
Loaded 5 loads with time series of shape (17568, 5)
Loaded 37 generators with time series of shape (17568, 18)
Loaded 9 lines
Loaded 8 links
Loaded 37 generators
Loaded 5 storage units
Total demand:  33504.55

Add carrier colours for plotting consistency.

Show the Python Code

# Add carrier colours for plotting consistency.
# Userved Energy bright red!

carrier_list = n.generators.carrier.unique()
carrier_colors = {
    "Biomass": '#127E2A',
    "Hydro": '#1E81D4',
    "Black Coal": "#39322D",
    "Solar": '#FDB324',
    "Rooftop Solar": '#FFE066',
    "Wind": '#3BBFE5',
    "Diesel": "#D486BA",
    "Brown Coal": "#715E50",
    "Gas": '#E6622D',
    "ROR": '#8ab2d4',
    "Battery": '#814ad4',
    "Pump Hydro": '#104775',
    "AC": '#999999',
    "DC": "#3277AF",
    "Unserved Energy": "#F40B16"
}

for carrier, color in carrier_colors.items():
    n.carriers.loc[carrier, 'color'] = color

Plot high-level NEM network (buses, lines & links).
Note: the lines/links below are not geographically accurate - just for representation.

Show the Python Code

# Plot high-level NEM network (buses, lines & links)

# Use PlateCarree projection for lat/lon
crs = ccrs.PlateCarree()

# Create figure and map
fig, ax = plt.subplots(figsize=(9, 6), subplot_kw={'projection': crs})

# Add base map features
ax.add_feature(cfeature.COASTLINE)
ax.add_feature(cfeature.BORDERS, linestyle=":")
ax.add_feature(cfeature.STATES, linewidth=0.5)
ax.set_extent([110, 155, -45, -10])  # Australia-wide view
# fig.patch.set_facecolor('#F0FFFF') 
# ax.set_facecolor('#F0FFFF')

# Plot buses
ax.scatter(
    n.buses["x"], n.buses["y"],
    color='red', s=50, transform=crs, zorder=5,
    label="Buses"
)

# Label each bus
for name, row in n.buses.iterrows():
    ax.text(
        row["x"], row["y"], name,
        transform=crs, fontsize=9, ha='left', va='bottom'
    )

# Plot lines (AC)
first_line = True
for _, line in n.lines.iterrows():
    x0, y0 = n.buses.loc[line.bus0, ["x", "y"]]
    x1, y1 = n.buses.loc[line.bus1, ["x", "y"]]
    color = carrier_colors.get(line.carrier, "gray")  # fallback if undefined
    
    ax.plot([x0, x1], [y0, y1], color=color, transform=crs, zorder=3,
            label="Lines (AC)" if first_line else None)
    first_line = False

# Plot links (DC)
first_link = True
for _, link in n.links.iterrows():
    x0, y0 = n.buses.loc[link.bus0, ["x", "y"]]
    x1, y1 = n.buses.loc[link.bus1, ["x", "y"]]
    color = carrier_colors.get(link.carrier, "blue")
    ax.plot([x0, x1], [y0, y1], color=color, linestyle="--", transform=crs, zorder=3,
            label="Links (DC)" if first_link else None)
    first_link = False

plt.title("PyPSA High-level NEM Network (Buses, Lines, Links)")
plt.legend()
plt.tight_layout()
plt.show()

03 Solve the network

Show the Python Code

# Solve the network
# Note: commented out to save run time..
# Change to open source solver if required, default is HiGHS:
# n.optimize()
#n.optimize(solver_name = "gurobi")

Show the Python Code

# Export the baseline network to NetCDF format so it is accessible for other tools.
# This is useful for sharing the network with others or for further analysis.
# n.export_to_netcdf("results/high-level_nem.nc")

Show the Python Code

# Read the network back from NetCDF format (if needed)
n = pypsa.Network("results/high-level_nem.nc")

INFO:pypsa.io:Imported network high-level_nem.nc has buses, carriers, generators, lines, links, loads, storage_units

.p_nom_opt is the optimised nominal capacity of a generator, after running n.optimize() in PyPSA. In the baseline model, it is expected that there will be zero unserved energy.

Show the Python Code

# Group by bus and carrier to sum the optimised nominal power in MW.
n.generators.groupby(["bus", "carrier"]).p_nom_opt.sum().reset_index() \
    .pivot(index='bus', columns='carrier', values='p_nom_opt') \
    .fillna(0).astype(int).sort_index()

carrier	Black Coal	Brown Coal	Diesel	Gas	Hydro	Rooftop Solar	Solar	Unserved Energy	Wind
bus
NSW1	8270	0	155	3278	7122	5494	7329	0	3212
QLD1	8119	0	453	3202	976	4604	5379	0	3402
SA1	0	0	454	3408	2	1885	735	0	2763
TAS1	0	0	48	372	2356	252	0	0	587
VIC1	0	4690	0	2491	2352	3664	2034	0	6044

04 Visualisations

Show the Python Code

# Plotting function for generation dispatch
# Interactive version (optional)


def plot_dispatch(n, time="2024", days=None, regions=None,
                   show_imports=True, show_curtailment=True,
                   scenario_name=None, scenario_objective=None, interactive=False):
     """
     Plot a generation dispatch stack by carrier for a PyPSA network, with optional
     net imports/exports and a region‑filtered curtailment overlay.

     Parameters
     ----------
     n : pypsa.Network
         The PyPSA network to plot.
     time : str, default "2024"
         Start of the time window (e.g. "2024", "2024-07", or "2024-07-15").
     days : int, optional
         Number of days from `time` to include in the plot.
     regions : list of str, optional
         Region bus names to filter by. If None, the entire network is included.
     show_imports : bool, default True
         Whether to include net imports/exports in the dispatch stack.
     show_curtailment : bool, default True
         Whether to calculate and plot VRE curtailment (Solar, Wind, Rooftop Solar).
     scenario_name : str, optional
         Scenario label to display below the title.
     scenario_objective : str, optional
         Objective description to display next to the legend.
     interactive : bool, default False
         Whether to create an interactive plot using Plotly instead of matplotlib.

     Notes
     -----
     - All power values are converted to GW.
     - Curtailment is plotted as a dashed black line if enabled.
     - Demand (load) is plotted as a solid green line.
     - Storage charging and net exports (negative values) are shown below zero.
     """
     
     
     if interactive:
         import plotly.graph_objects as go
         from plotly.subplots import make_subplots
         import plotly.express as px
     else:
         import matplotlib.pyplot as plt
     
     # 1) REGION MASKS
     if regions is not None:
         gen_mask   = n.generators.bus.isin(regions)
         sto_mask   = n.storage_units.bus.isin(regions) if not n.storage_units.empty else []
         store_mask = n.stores.bus.isin(regions) if not n.stores.empty else []
         region_buses = set(regions)
         
     else:
        gen_mask = pd.Series(True, index=n.generators.index)
        sto_mask = pd.Series(True, index=n.storage_units.index) if not n.storage_units.empty else pd.Series(dtype=bool)
        store_mask = pd.Series(True, index=n.stores.index) if not n.stores.empty else pd.Series(dtype=bool)
        region_buses = set(n.buses.index)

     # 2) AGGREGATE BY CARRIER (GW)
     def _agg(df_t, df_stat, mask):
         return (
             df_t.loc[:, mask]
                 .T
                 .groupby(df_stat.loc[mask, 'carrier'])
                 .sum()
                 .T
                 .div(1e3)
         )

     p_by_carrier = _agg(n.generators_t.p, n.generators, gen_mask)
     if not n.storage_units.empty:
         p_by_carrier = pd.concat([p_by_carrier,
                                   _agg(n.storage_units_t.p, n.storage_units, sto_mask)],
                                  axis=1)
     if not n.stores.empty:
         p_by_carrier = pd.concat([p_by_carrier,
                                   _agg(n.stores_t.p, n.stores, store_mask)],
                                  axis=1)

     # 3) TIME WINDOW
     parts = time.split("-")
     if len(parts) == 1:
         start = pd.to_datetime(f"{parts[0]}-01-01")
     elif len(parts) == 2:
         start = pd.to_datetime(f"{parts[0]}-{parts[1]}-01")
     else:
         start = pd.to_datetime(time)

     if days is not None:
         end = start + pd.Timedelta(days=days) - pd.Timedelta(hours=1)
     elif len(parts) == 1:
         end = pd.to_datetime(f"{parts[0]}-12-31 23:00")
     elif len(parts) == 2:
         end = start + pd.offsets.MonthEnd(0) + pd.Timedelta(hours=23)
     else:
         end = start + pd.Timedelta(hours=23)

     p_slice = p_by_carrier.loc[start:end].copy()
     # drop carriers with zero activity
     zero = p_slice.columns[p_slice.abs().sum() == 0]
     p_slice.drop(columns=zero, inplace=True)

     # 4) IMPORTS/EXPORTS
     if show_imports:
         ac = ( n.lines_t.p0.loc[start:end, n.lines.bus1.isin(region_buses) & ~n.lines.bus0.isin(region_buses)].sum(axis=1)
              + n.lines_t.p1.loc[start:end, n.lines.bus0.isin(region_buses) & ~n.lines.bus1.isin(region_buses)].sum(axis=1) )
         dc = ( n.links_t.p0.loc[start:end, n.links.bus1.isin(region_buses) & ~n.links.bus0.isin(region_buses)].sum(axis=1)
              + n.links_t.p1.loc[start:end, n.links.bus0.isin(region_buses) & ~n.links.bus1.isin(region_buses)].sum(axis=1) )
         p_slice['Imports/Exports'] = (ac + dc).div(1e3)
         if 'Imports/Exports' not in n.carriers.index:
             n.carriers.loc['Imports/Exports','color']='#7f7f7f'

     # 5) LOAD SERIES
     if regions:
         load_cols = [c for c in n.loads[n.loads.bus.isin(regions)].index if c in n.loads_t.p_set]
         load_series = n.loads_t.p_set[load_cols].sum(axis=1)
     else:
         load_series = n.loads_t.p_set.sum(axis=1)
     load_series = load_series.loc[start:end].div(1e3)

     # 6) VRE CURTAILMENT (GW) if requested
     if show_curtailment:
         vre = ['Solar','Wind', 'Rooftop Solar']
         mask_vre = gen_mask & n.generators.carrier.isin(vre)
         avail = (n.generators_t.p_max_pu.loc[start:end, mask_vre]
                  .multiply(n.generators.loc[mask_vre,'p_nom'], axis=1))
         disp  = n.generators_t.p.loc[start:end, mask_vre]
         curtail = (avail.sub(disp, fill_value=0)
                        .clip(lower=0)
                        .sum(axis=1)
                        .div(1e3))
     else:
         curtail = None

     # 7) PLOT
     title_tail = f" for {', '.join(regions)}" if regions else ''
     plot_title = f"Dispatch by Carrier: {start.date()} to {end.date()}{title_tail}"
     
     if interactive:
         # PLOTLY INTERACTIVE PLOT
         fig = go.Figure()
         
         fig.update_layout(
            plot_bgcolor='#F0FFFF',
            xaxis=dict(gridcolor='#DDDDDD'),
            yaxis=dict(gridcolor='#DDDDDD')        # Plot area background
         ) 

         # Prepare data for stacked area plot
         positive_data = p_slice.where(p_slice > 0).fillna(0)
         negative_data = p_slice.where(p_slice < 0).fillna(0)
         
         # Add positive generation as stacked area
         for i, col in enumerate(positive_data.columns):
             if positive_data[col].sum() > 0:
                 color = n.carriers.loc[col, 'color']
                 # Only add points where value > 0.001
                 mask = positive_data[col].abs() > 0.001
                 if mask.any():
                     fig.add_trace(go.Scatter(
                         x=positive_data.index[mask],
                         y=positive_data[col][mask],
                         mode='lines',
                         fill='tonexty' if i > 0 else 'tozeroy',
                         line=dict(width=0, color=color),
                         fillcolor=color,
                         name=col,
                         stackgroup='positive',
                         hovertemplate='<b>%{fullData.name}</b><br>Power: %{y:.3f} GW<extra></extra>',
                         showlegend=True
                     ))
         
         # Add negative generation (storage charging, exports)
         for col in negative_data.columns:
             if negative_data[col].sum() < 0:
                 color = n.carriers.loc[col, 'color']
                 # Only add points where value < -0.001
                 mask = negative_data[col].abs() > 0.001
                 if mask.any():
                     fig.add_trace(go.Scatter(
                         x=negative_data.index[mask],
                         y=negative_data[col][mask],
                         mode='lines',
                         fill='tonexty',
                         line=dict(width=0, color=color),
                         fillcolor=color,
                         name=col,
                         stackgroup='negative',
                         hovertemplate='<b>%{fullData.name}</b><br>Power: %{y:.2f} GW<extra></extra>',
                         showlegend=True
                     ))
         
         # Add demand line (always show)
         fig.add_trace(go.Scatter(
             x=load_series.index,
             y=load_series,
             mode='lines',
             line=dict(color='green', width=2),
             name='Demand',
             hovertemplate='<b>Demand</b><br>Power: %{y:.2f} GW<extra></extra>',
             showlegend=True
         ))
         
         # Add curtailment line if requested
         if show_curtailment and curtail is not None:
             fig.add_trace(go.Scatter(
                 x=curtail.index,
                 y=curtail,
                 mode='lines',
                 line=dict(color='black', width=2, dash='dash'),
                 name='Curtailment',
                 hovertemplate='<b>Curtailment</b><br>Power: %{y:.2f} GW<extra></extra>',
                 showlegend=True
                 ))
         
         # Update layout
         fig.update_layout(
             title=plot_title,
             xaxis_title='Time',
             yaxis_title='Power (GW)',
             hovermode='x unified',
             hoverlabel=dict(
                 bgcolor="white",
                 bordercolor="black",
                 font_size=12,
             ),
             legend=dict(
                 x=1.02,
                 y=1,
                 bgcolor='rgba(255,255,255,0.8)',
                 bordercolor='rgba(0,0,0,0.2)',
                 borderwidth=1
             ),
             width=780,
             height=500
         )
         
         # Add scenario annotations
         annotations = []
         if scenario_name:
             annotations.append(
                 dict(
                     x=1.02, y=-0.05,
                     xref='paper', yref='paper',
                     text=f"Scenario: {scenario_name}",
                     showarrow=False,
                     font=dict(size=10, color='gray'),
                     xanchor='center',
                     yanchor='top'
                 )
             )
         
         if annotations:
             fig.update_layout(annotations=annotations)
         
         fig.show()
        #  return fig
         
     else:
         # MATPLOTLIB STATIC PLOT 
         fig, ax = plt.subplots(figsize=(8.4, 6.5)) #12,6.5
         cols = p_slice.columns.map(lambda c: n.carriers.loc[c,'color'])
         p_slice.where(p_slice>0).plot.area(ax=ax,linewidth=0,color=cols)
         neg = p_slice.where(p_slice<0).dropna(how='all',axis=1)
         if not neg.empty:
             neg_cols=[n.carriers.loc[c,'color'] for c in neg.columns]
             neg.plot.area(ax=ax,linewidth=0,color=neg_cols)
         load_series.plot(ax=ax,color='g',linewidth=1.5,label='Demand')
         if show_curtailment and curtail is not None:
             curtail.plot(ax=ax,color='k',linestyle='--',linewidth=1.2,label='Curtailment')

         # limits & legend
         up = max(p_slice.where(p_slice>0).sum(axis=1).max(),
                  load_series.max(),
                  curtail.max() if curtail is not None else 0)
         dn = min(p_slice.where(p_slice<0).sum(axis=1).min(), load_series.min())
         ax.set_ylim(dn if not np.isclose(up,dn) else dn-0.1, up)
        #  fig.patch.set_facecolor('#F0FFFF') 
         ax.set_facecolor('#F0FFFF')
         h,l = ax.get_legend_handles_labels()
         seen={} ; fh,fl=[],[]
         for hh,ll in zip(h,l):
             if ll not in seen: fh.append(hh);fl.append(ll);seen[ll]=True
         ax.legend(fh,fl,loc=(1.02,0.67), fontsize=9)

         # scenario text
         if scenario_objective:
             ax.text(1.02,0.01,f"Objective:\n{scenario_objective}",transform=ax.transAxes,
                     fontsize=8,va='bottom',ha='left',bbox=dict(facecolor='white',alpha=0.7,edgecolor='none'))
         if scenario_name:
             ax.text(1.02,-0.05,f"Scenario: {scenario_name}",transform=ax.transAxes,
                     fontsize=9,color='gray',ha='center',va='top')

         ax.set_ylabel('GW')
         ax.set_title(plot_title)
         plt.tight_layout()
         plt.show()

`0_2024_baseline` scenario

Review dispatch for a 5-day period in NSW.

Note here that very little gas is dispatched due to having a higher marginal price than black coal. In reality the ramping up and down of coal plants is smoother than this. The 30-min resolution does not constrain the ramping much compared to actual five-minute dispatching in the NEM. For example, Eraring black coal units maximum ramp-up of 4MW/min, is 120MW per model snapshot³.

Show the Python Code

# plot_dispatch function example - static plot of NSW1 region for 5 days starting from 2024-07-01
plot_dispatch(n, time="2024-07-01", days=5, regions=["NSW1"], show_imports=True)

The interactive version is good for exploratory data analysis.

Show the Python Code

# plot_dispatch function example - interactive plot of VIC1 region for 5 days starting from 2024-01-19
plot_dispatch(n, time="2024-01-10", days=5, regions=["VIC1"], show_imports=True, interactive=True)

Review the 2024 baseline model generator capacities.

Show the Python Code

def plot_generator_capacity_by_carrier(network):
    """
    Plot total generator capacity by carrier for a given PyPSA network,
    excluding carriers with zero total capacity, and colour bars by carrier colour.
    """
    # 1) sum and filter
    capacity_by_carrier = (
        network.generators
               .groupby("carrier")
               .p_nom_opt
               .sum()
               .div(1e3)
    )
    capacity_by_carrier = capacity_by_carrier[capacity_by_carrier > 0]

    # 2) get colors
    colors = [
        network.carriers.loc[c, 'color']
            if c in network.carriers.index and 'color' in network.carriers.columns
            else 'gray'
        for c in capacity_by_carrier.index
    ]

    # 3) create fig & ax, set backgrounds
    fig, ax = plt.subplots(figsize=(8, 5))
    # fig.patch.set_facecolor('#F0FFFF')    # full-figure bg
               # axes bg

    # 4) plot onto that ax
    capacity_by_carrier.plot.barh(ax=ax, color=colors)

    # 5) labels & layout
    ax.set_facecolor('#F0FFFF')
    ax.set_xlabel("GW")
    ax.set_title("Total Generator Capacity by Carrier")
    fig.tight_layout()
    plt.show()

plot_generator_capacity_by_carrier(n)

According to the AEMO’s ISP Step Change - all coal fired generation will be retired by 2038. This is ‘just around the corner’ in terms of long lead times for generators and associated infrastructure build-outs.

Plot peak demand vs installed generation capacity per region.

Show the Python Code

def plot_peak_demand_vs_capacity(network, stack_carrier=False, relative=False, sort_by="Max Demand"):
    """
    Plot peak demand vs installed generation capacity per bus, filtering out zero-capacity carriers.
    Consistent styling and behavior with plot_total_demand_vs_generation function.

    Parameters:
    -----------
    network : pypsa.Network
        The PyPSA network object.
    stack_carrier : bool, default=False
        Whether to stack capacity by carrier or plot total capacity.
    relative : bool, default=False
        If True, plot capacity as a percentage of peak demand.
    sort_by : str, optional
        Column to sort buses by on the x-axis. Options: "Max Demand", "Total Capacity", or None.
    """
    # Define preferred carrier order (thermal/dispatchable at bottom, renewables on top)
    preferred_order = [
        'Black Coal', 'Brown Coal', 'Gas', 'Diesel', 'Hydro', 
        'Wind', 'Solar', 'Rooftop Solar', 'Unserved Energy'
    ]
    
    # --- 1) Max demand per bus ---
    max_demand = network.loads_t.p_set.max()

    # --- 2) Capacity per bus and carrier ---
    capacity_by_bus_carrier = (
        network.generators.groupby(['bus', 'carrier'])
        .p_nom_opt.sum()
        .unstack(fill_value=0)
    )
    # Filter out carriers with zero total capacity
    nonzero_carriers = capacity_by_bus_carrier.columns[capacity_by_bus_carrier.sum(axis=0) > 0]
    capacity_by_bus_carrier = capacity_by_bus_carrier[nonzero_carriers]
    total_capacity = capacity_by_bus_carrier.sum(axis=1)

    # --- 3) Combine DataFrame and filter out zero-demand & zero-capacity buses ---
    df = pd.concat([max_demand.rename("Max Demand"),
                    total_capacity.rename("Total Capacity")], axis=1).fillna(0)
    df = df[(df["Max Demand"] > 0) | (df["Total Capacity"] > 0)]

    # --- 4) Relative scaling if requested ---
    if relative:
        # Create a copy to avoid modifying original data
        df_plot = df.copy()
        # avoid div-by-zero
        df_plot["Max Demand"] = df_plot["Max Demand"].replace(0, np.nan)
        relative_capacity = capacity_by_bus_carrier.div(df["Max Demand"], axis=0) * 100
        relative_capacity = relative_capacity.fillna(0)  # Handle NaN from division by zero
        df_plot["Total Capacity"] = df_plot["Total Capacity"] / df_plot["Max Demand"] * 100
        df_plot["Max Demand"] = 100
        ylabel = "Capacity as % of Peak Demand"
    else:
        df_plot = df.copy()
        # convert to GW for consistency
        df_plot[["Max Demand", "Total Capacity"]] = df_plot[["Max Demand", "Total Capacity"]] / 1e3
        relative_capacity = capacity_by_bus_carrier / 1e3
        ylabel = "Power (GW)"

    # --- 5) Sort consistently (always descending, handle NaN) ---
    if sort_by in df_plot.columns:
        df_plot = df_plot.sort_values(by=sort_by, ascending=False, na_position='last')
        
    # Reindex capacity data to match sorted order
    relative_capacity = relative_capacity.reindex(df_plot.index, fill_value=0)

    # --- 6) Plotting with consistent styling ---
    fig, ax = plt.subplots(figsize=(8, 5))  # Consistent figure size
    bar_width = 0.35
    bar_pos = np.arange(len(df_plot))

    if stack_carrier:
        # Stack each non-zero carrier in preferred order
        bottom = np.zeros(len(df_plot))
        
        # Get available carriers and order them according to preferred_order
        available_carriers = list(capacity_by_bus_carrier.columns)
        
        # Order carriers: first by preferred order, then any remaining carriers
        ordered_carriers = []
        for carrier in preferred_order:
            if carrier in available_carriers:
                ordered_carriers.append(carrier)
        
        # Add any carriers not in preferred_order at the end
        for carrier in available_carriers:
            if carrier not in ordered_carriers:
                ordered_carriers.append(carrier)
        
        carriers = ordered_carriers
        
        for carrier in carriers:
            vals = relative_capacity[carrier].fillna(0).values
            
            # Get color from network.carriers if available, otherwise use default colors
            if (hasattr(network, 'carriers') and 
                carrier in network.carriers.index and 
                'color' in network.carriers.columns):
                color = network.carriers.loc[carrier, 'color']
            else:
                # Default color cycle
                colors = plt.cm.Set3(np.linspace(0, 1, len(carriers)))
                color = colors[carriers.index(carrier)]
                
            ax.bar(bar_pos + bar_width/2, vals, bar_width,
                   bottom=bottom, label=f'Capacity ({carrier})', color=color)
            bottom += vals
            
        # Plot peak demand on left
        ax.bar(bar_pos - bar_width/2, df_plot["Max Demand"], bar_width,
               label='Peak Demand', color='gray', alpha=0.7)
    else:
        ax.bar(bar_pos - bar_width/2, df_plot["Max Demand"], bar_width,
               label='Peak Demand', color='gray', alpha=0.7)
        ax.bar(bar_pos + bar_width/2, df_plot["Total Capacity"], bar_width,
               label='Total Capacity', color='tab:blue')

    # --- 7) Consistent labels and styling ---
    ax.set_xticks(bar_pos)
    ax.set_xticklabels(df_plot.index, rotation=45, ha='right')  # Consistent rotation
    ax.set_ylabel(ylabel)
    ax.set_xlabel("Region")
    ax.set_title("Peak Demand vs Generation Capacity per Region" + (" (Relative)" if relative else ""))
    ax.set_facecolor('#F0FFFF')  # Consistent background color
    ax.legend(loc='upper right')  # Consistent legend placement
    ax.grid(True, alpha=0.3)  # Consistent grid styling
    plt.tight_layout()
    plt.show()
    
    # Return the actual values (not plot data which might be relative)
    # return df


# Call the function with your preferred ordering
plot_peak_demand_vs_capacity(n, stack_carrier=True)

Plot total electricity demand vs generation per bus/region.

Show the Python Code

def plot_total_demand_vs_generation(network, stack_carrier=False, relative=False):
    """
    Plot total electricity demand vs generation per bus.
    Properly accounts for snapshot weightings to handle different temporal resolutions.

    Parameters:
    - network: PyPSA Network object
    - stack_carrier: If True, stack generation by carrier (color-coded)
    - relative: If True, show generation as % of total demand (demand = 100%)
    """
    
    # Define preferred carrier order (thermal/dispatchable at bottom, renewables on top)
    preferred_order = [
        'Black Coal', 'Brown Coal', 'Gas', 'Diesel', 'Hydro', 
        'Wind', 'Solar', 'Rooftop Solar', 'Unserved Energy'
    ]
    
    def get_snapshot_weights(network):
        """Get snapshot weightings in a robust way"""
        if hasattr(network.snapshot_weightings, 'generators'):
            return network.snapshot_weightings['generators']
        elif 'generators' in network.snapshot_weightings.columns:
            return network.snapshot_weightings['generators']
        else:
            # If single column or series, use first column
            return network.snapshot_weightings.iloc[:, 0] if network.snapshot_weightings.ndim > 1 else network.snapshot_weightings
    
    # Get snapshot weights
    weights = get_snapshot_weights(network)
    
    # Total demand per bus in GWh (properly weighted)
    demand_energy = network.loads_t.p_set.multiply(weights, axis=0)
    total_demand_per_bus = demand_energy.sum().div(1e3)
    total_demand_per_bus.name = "Total Demand"

    # Total generation per generator in GWh (properly weighted)
    gen_energy_weighted = network.generators_t.p.multiply(weights, axis=0)
    gen_energy = gen_energy_weighted.sum().div(1e3)
    
    gen_info = network.generators[["bus", "carrier"]]
    gen_energy_by_carrier = (
        gen_info.assign(energy=gen_energy)
        .groupby(["bus", "carrier"])["energy"]
        .sum()
        .unstack(fill_value=0)
    )
    total_generation_per_bus = gen_energy_by_carrier.sum(axis=1)
    total_generation_per_bus.name = "Total Generation"

    # Join and filter
    generation_vs_demand = pd.concat([total_demand_per_bus, total_generation_per_bus], axis=1).fillna(0)
    generation_vs_demand = generation_vs_demand.loc[
        (generation_vs_demand["Total Demand"] > 0) | (generation_vs_demand["Total Generation"] > 0)
    ]

    if relative:
        # Create a copy to avoid modifying original data
        generation_vs_demand_plot = generation_vs_demand.copy()
        generation_vs_demand_plot["Total Demand"].replace(0, np.nan, inplace=True)  # avoid div by 0
        generation_vs_demand_plot = generation_vs_demand_plot.sort_values(by="Total Demand", ascending=False, na_position='last')
        
        # Calculate relative generation by carrier
        relative_generation = gen_energy_by_carrier.div(generation_vs_demand["Total Demand"], axis=0) * 100
        relative_generation = relative_generation.fillna(0)  # Handle NaN from division by zero
        
        generation_vs_demand_plot["Total Generation"] = (
            generation_vs_demand_plot["Total Generation"] / generation_vs_demand_plot["Total Demand"] * 100
        )
        generation_vs_demand_plot["Total Demand"] = 100
        ylabel = "Generation vs Demand (%)"
    else:
        generation_vs_demand_plot = generation_vs_demand.copy()
        relative_generation = gen_energy_by_carrier
        ylabel = "Energy (GWh)"

    # Sort by Total Demand (descending) and reindex    
    generation_vs_demand_plot = generation_vs_demand_plot.sort_values(by="Total Demand", ascending=False, na_position='last')
    relative_generation = relative_generation.reindex(generation_vs_demand_plot.index, fill_value=0)

    # Plot
    fig, ax = plt.subplots(figsize=(8, 5))
    bar_width = 0.35
    bar_positions = np.arange(len(generation_vs_demand_plot))

    if stack_carrier:
        bottom = np.zeros(len(generation_vs_demand_plot))
        
        # Get carriers that actually have generation
        available_carriers = [
            c for c in relative_generation.columns
            if (relative_generation[c].sum() > 0)
        ]
        
        # Order carriers according to preferred_order
        ordered_carriers = []
        for carrier in preferred_order:
            if carrier in available_carriers:
                ordered_carriers.append(carrier)
        
        # Add any carriers not in preferred_order at the end
        for carrier in available_carriers:
            if carrier not in ordered_carriers:
                ordered_carriers.append(carrier)
        
        carriers = ordered_carriers
        
        for carrier in carriers:
            values = relative_generation[carrier].fillna(0)
            
            # Get color from network.carriers if available, otherwise use default colors
            if (hasattr(network, 'carriers') and 
                carrier in network.carriers.index and 
                'color' in network.carriers.columns):
                color = network.carriers.loc[carrier, 'color']
            else:
                # Default color cycle
                colors = plt.cm.Set3(np.linspace(0, 1, len(carriers)))
                color = colors[carriers.index(carrier)]
            
            ax.bar(
                bar_positions + bar_width / 2, values, bar_width,
                label=f'Generation ({carrier})', bottom=bottom, color=color
            )
            bottom += values.values

        ax.bar(
            bar_positions - bar_width / 2,
            generation_vs_demand_plot["Total Demand"],
            bar_width,
            label="Total Demand",
            color="gray", alpha=0.7
        )
    else:
        ax.bar(
            bar_positions - bar_width / 2,
            generation_vs_demand_plot["Total Demand"],
            bar_width,
            label="Total Demand",
            color="gray", alpha=0.7
        )
        ax.bar(
            bar_positions + bar_width / 2,
            generation_vs_demand_plot["Total Generation"],
            bar_width,
            label="Total Generation",
            color="tab:blue"
        )

    ax.set_xlabel("Region")
    ax.set_ylabel(ylabel)
    ax.set_title("Total Demand vs Total Generation per Region" + (" (Relative)" if relative else ""))
    ax.set_xticks(bar_positions)
    ax.set_xticklabels(generation_vs_demand_plot.index, rotation=45, ha="right")
    ax.set_facecolor('#F0FFFF')
    ax.legend(loc='upper right')
    ax.grid(True, alpha=0.3)
    plt.tight_layout()
    plt.show()
    
    # Return the actual energy values (not the plot data which might be relative)
    # return generation_vs_demand

# Call the function with your preferred ordering
plot_total_demand_vs_generation(n, stack_carrier=True, relative=False)

Load previously generated scenario data.

The 0_2024_baseline scenario is the baseline and hence scale-up factor is 1. The 19 others are scaled up from the baseline model according to the Objective column.

Note: some of the scenarios have unserved energy (bold, red font), which tends to show up in NSW and VIC before other regions.

Show the Python Code

# Previously saved results generated from the generate_scenarios function in the Appendix
# Note: 0_2024_baseline is the baseline scenario and is equivalent to 'high-level_nem'.
df_results = pd.read_csv("results/scenarios/scenarios_summary_20251013_1439.csv")

# Format in a more user-friendly format for HTML
results_tbl = df_results.copy()
for col in df_results.select_dtypes(include=['object']).columns:
    results_tbl[col] = results_tbl[col].str.replace('\n', '<br>', regex=False)

# Drop to make more space.
# 'Wind & Solar Curtailment (GWh)'
results_tbl = results_tbl.drop(['Unserved Energy (GWh)'], axis=1)

# Get the last column name
last_col = results_tbl.columns[-1]

# Apply conditional formatting to the last column (which contains dictionary-like data)
def format_last_column(val):
    if pd.isna(val) or val == '':
        return str(val)
    
    import re
    
    # Use regex to find all key: value patterns
    pattern = r'(\w+\d*): ([\d.-]+)'
    
    def replace_positive_values(match):
        key = match.group(1)
        value = match.group(2)
        try:
            num_val = float(value)
            if num_val > 0:
                return f'{key}: <span style="color: red; font-weight: bold;">{value}</span>'
            else:
                return f'{key}: {value}'
        except (ValueError, TypeError):
            return f'{key}: {value}'
    
    # Replace all matches with formatted versions
    result = re.sub(pattern, replace_positive_values, str(val))
    return result

# Apply formatting to the last column
results_tbl[last_col] = results_tbl[last_col].apply(format_last_column)

html_table = results_tbl.to_html(escape=False, index=False)

# Add container div with fixed height and scrolling
html_table = f'''
<div style="height: 600px; overflow-y: auto; border: 1px solid #ddd;">
{html_table}
</div>
'''

# Style the table with sticky header
html_table = html_table.replace('<table', 
    '<table style="width: 100%; table-layout: fixed; font-size: 11px; font-family: monospace; border-collapse: collapse;"')
html_table = html_table.replace('<th', 
    '<th style="position: sticky; top: 0; background-color: #f5f5f5; border: 1px solid #ddd; padding: 8px; z-index: 10;"')
html_table = html_table.replace('<td', 
    '<td style="min-width: 250px; word-wrap: break-word; padding: 4px; border: 1px solid #eee;"')

print(html_table)


Scenario	Objective	GPG (GWh)	Battery Capacity (GW)	Generator Capacity (GW)	Wind & Solar Curtailment (GWh)	Unserved by Region (GWh)
0_2024_baseline	NSW1-Black Coal: x1.0 NSW1-Rooftop Solar: x1.0 NSW1-Solar: x1.0 NSW1-Wind: x1.0 NSW1-Battery: x1.0 QLD1-Black Coal: x1.0 QLD1-Rooftop Solar: x1.0 QLD1-Solar: x1.0 QLD1-Wind: x1.0 QLD1-Battery: x1.0 SA1-Rooftop Solar: x1.0 SA1-Solar: x1.0 SA1-Wind: x1.0 SA1-Battery: x1.0 VIC1-Brown Coal: x1.0 VIC1-Rooftop Solar: x1.0 VIC1-Solar: x1.0 VIC1-Wind: x1.0 VIC1-Battery: x1.0	1663.143	12.101	Black Coal: 16389.0 Brown Coal: 4690.0 Diesel: 1110.0 Gas: 12751.0 Hydro: 12808.0 Rooftop Solar: 15901.42 Solar: 15477.0 Wind: 16008.0	4607.511	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
1_BalancedTransition	NSW1-Black Coal: x0.5 NSW1-Rooftop Solar: x1.2 NSW1-Solar: x1.2 NSW1-Wind: x1.25 NSW1-Battery: x1.5 QLD1-Black Coal: x0.5 QLD1-Rooftop Solar: x1.2 QLD1-Solar: x1.2 QLD1-Wind: x1.25 QLD1-Battery: x1.5 SA1-Rooftop Solar: x1.2 SA1-Solar: x1.2 SA1-Wind: x1.25 SA1-Battery: x1.5 TAS1-Wind: x1.25 TAS1-Battery: x1.5 VIC1-Brown Coal: x0.5 VIC1-Rooftop Solar: x1.2 VIC1-Solar: x1.2 VIC1-Wind: x1.25 VIC1-Battery: x1.5	6290.599	18.152	Black Coal: 8194.5 Brown Coal: 2345.0 Diesel: 1110.0 Gas: 12751.0 Hydro: 12808.0 Rooftop Solar: 19031.214 Solar: 18572.4 Wind: 20010.0	7362.259	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
2_BalancedAggressiveTransition	NSW1-Black Coal: x0.25 NSW1-Rooftop Solar: x1.25 NSW1-Solar: x1.25 NSW1-Wind: x1.3 NSW1-Battery: x2.0 QLD1-Black Coal: x0.25 QLD1-Rooftop Solar: x1.25 QLD1-Solar: x1.25 QLD1-Wind: x1.3 QLD1-Battery: x2.0 SA1-Rooftop Solar: x1.25 SA1-Solar: x1.25 SA1-Wind: x1.3 SA1-Battery: x2.0 TAS1-Wind: x1.3 TAS1-Battery: x2.0 VIC1-Brown Coal: x0.25 VIC1-Rooftop Solar: x1.25 VIC1-Solar: x1.25 VIC1-Wind: x1.3 VIC1-Battery: x2.0	21090.899	24.202	Black Coal: 4097.25 Brown Coal: 1172.5 Diesel: 1110.0 Gas: 12751.0 Hydro: 12808.0 Rooftop Solar: 19813.662 Solar: 19346.25 Wind: 20810.4	5551.653	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
3.0_VreStorageRampTransition	NSW1-Black Coal: x0.15 NSW1-Rooftop Solar: x1.3 NSW1-Solar: x1.3 NSW1-Wind: x1.5 NSW1-Battery: x4.0 QLD1-Black Coal: x0.15 QLD1-Rooftop Solar: x1.3 QLD1-Solar: x1.3 QLD1-Wind: x1.5 QLD1-Battery: x4.0 SA1-Rooftop Solar: x1.3 SA1-Solar: x1.3 SA1-Wind: x1.5 SA1-Battery: x4.0 TAS1-Wind: x1.5 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.15 VIC1-Rooftop Solar: x1.3 VIC1-Solar: x1.3 VIC1-Wind: x1.5 VIC1-Battery: x4.0	19321.371	48.404	Black Coal: 2458.35 Brown Coal: 703.5 Diesel: 1110.0 Gas: 12751.0 Hydro: 12808.0 Rooftop Solar: 20596.111 Solar: 20120.1 Wind: 24012.0	1340.579	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
3.1_VreStorageRampGasReduction	NSW1-Black Coal: x0.15 NSW1-Gas: x0.9 NSW1-Rooftop Solar: x1.3 NSW1-Solar: x1.3 NSW1-Wind: x1.5 NSW1-Battery: x4.0 QLD1-Black Coal: x0.15 QLD1-Gas: x0.9 QLD1-Rooftop Solar: x1.3 QLD1-Solar: x1.3 QLD1-Wind: x1.5 QLD1-Battery: x4.0 SA1-Gas: x0.9 SA1-Rooftop Solar: x1.3 SA1-Solar: x1.3 SA1-Wind: x1.5 SA1-Battery: x4.0 TAS1-Gas: x0.9 TAS1-Wind: x1.5 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.15 VIC1-Gas: x0.9 VIC1-Rooftop Solar: x1.3 VIC1-Solar: x1.3 VIC1-Wind: x1.5 VIC1-Battery: x4.0	19247.624	48.404	Black Coal: 2458.35 Brown Coal: 703.5 Diesel: 1110.0 Gas: 11475.9 Hydro: 12808.0 Rooftop Solar: 20596.111 Solar: 20120.1 Wind: 24012.0	1330.714	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
3.1.1_VreStorageRampGasReduction	NSW1-Black Coal: x0.15 NSW1-Gas: x0.9 NSW1-Rooftop Solar: x1.3 NSW1-Solar: x1.3 NSW1-Wind: x1.5 NSW1-Battery: x5.0 QLD1-Black Coal: x0.15 QLD1-Gas: x0.9 QLD1-Rooftop Solar: x1.3 QLD1-Solar: x1.3 QLD1-Wind: x1.5 QLD1-Battery: x4.0 SA1-Gas: x0.9 SA1-Rooftop Solar: x1.3 SA1-Solar: x1.3 SA1-Wind: x1.5 SA1-Battery: x4.0 TAS1-Gas: x0.9 TAS1-Wind: x1.5 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.15 VIC1-Gas: x0.9 VIC1-Rooftop Solar: x1.3 VIC1-Solar: x1.3 VIC1-Wind: x1.5 VIC1-Battery: x5.0	18378.667	55.103	Black Coal: 2458.35 Brown Coal: 703.5 Diesel: 1110.0 Gas: 11475.9 Hydro: 12808.0 Rooftop Solar: 20596.111 Solar: 20120.1 Wind: 24012.0	1150.139	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
3.1.2_VreStorageRampGasReduction	NSW1-Black Coal: x0.15 NSW1-Gas: x0.9 NSW1-Rooftop Solar: x1.3 NSW1-Solar: x1.3 NSW1-Wind: x1.5 NSW1-Battery: x6.0 QLD1-Black Coal: x0.15 QLD1-Gas: x0.9 QLD1-Rooftop Solar: x1.3 QLD1-Solar: x1.3 QLD1-Wind: x1.5 QLD1-Battery: x4.0 SA1-Gas: x0.9 SA1-Rooftop Solar: x1.3 SA1-Solar: x1.3 SA1-Wind: x1.5 SA1-Battery: x4.0 TAS1-Gas: x0.9 TAS1-Wind: x1.5 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.15 VIC1-Gas: x0.9 VIC1-Rooftop Solar: x1.3 VIC1-Solar: x1.3 VIC1-Wind: x1.5 VIC1-Battery: x6.0	17921.405	61.802	Black Coal: 2458.35 Brown Coal: 703.5 Diesel: 1110.0 Gas: 11475.9 Hydro: 12808.0 Rooftop Solar: 20596.111 Solar: 20120.1 Wind: 24012.0	1057.004	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
3.2_VreStorageRampGasReduction	NSW1-Black Coal: x0.15 NSW1-Gas: x0.75 NSW1-Rooftop Solar: x1.3 NSW1-Solar: x1.3 NSW1-Wind: x1.5 NSW1-Battery: x4.0 QLD1-Black Coal: x0.15 QLD1-Gas: x0.75 QLD1-Rooftop Solar: x1.3 QLD1-Solar: x1.3 QLD1-Wind: x1.5 QLD1-Battery: x4.0 SA1-Gas: x0.75 SA1-Rooftop Solar: x1.3 SA1-Solar: x1.3 SA1-Wind: x1.5 SA1-Battery: x4.0 TAS1-Gas: x0.75 TAS1-Wind: x1.5 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.15 VIC1-Gas: x0.75 VIC1-Rooftop Solar: x1.3 VIC1-Solar: x1.3 VIC1-Wind: x1.5 VIC1-Battery: x4.0	18902.276	48.404	Black Coal: 2458.35 Brown Coal: 703.5 Diesel: 1110.0 Gas: 9563.25 Hydro: 12808.0 Rooftop Solar: 20596.111 Solar: 20120.1 Wind: 24012.0	1334.508	NSW1: 52.418 QLD1: 6.951 SA1: 0.0 TAS1: 0.875 VIC1: 211.738
3.3_VreStorageRampGasReduction	NSW1-Black Coal: x0.15 NSW1-Gas: x0.5 NSW1-Rooftop Solar: x1.3 NSW1-Solar: x1.3 NSW1-Wind: x1.5 NSW1-Battery: x4.0 QLD1-Black Coal: x0.15 QLD1-Gas: x0.5 QLD1-Rooftop Solar: x1.3 QLD1-Solar: x1.3 QLD1-Wind: x1.5 QLD1-Battery: x4.0 SA1-Gas: x0.5 SA1-Rooftop Solar: x1.3 SA1-Solar: x1.3 SA1-Wind: x1.5 SA1-Battery: x4.0 TAS1-Gas: x0.5 TAS1-Wind: x1.5 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.15 VIC1-Gas: x0.5 VIC1-Rooftop Solar: x1.3 VIC1-Solar: x1.3 VIC1-Wind: x1.5 VIC1-Battery: x4.0	16700.420	48.404	Black Coal: 2458.35 Brown Coal: 703.5 Diesel: 1110.0 Gas: 6375.5 Hydro: 12808.0 Rooftop Solar: 20596.111 Solar: 20120.1 Wind: 24012.0	1298.623	NSW1: 1689.721 QLD1: 113.523 SA1: 9.647 TAS1: 27.268 VIC1: 1084.068
4_4xVreTransitionZeroCoal	NSW1-Black Coal: x0.0 NSW1-Rooftop Solar: x2.0 NSW1-Solar: x4.0 NSW1-Wind: x4.0 NSW1-Battery: x4.0 QLD1-Black Coal: x0.0 QLD1-Rooftop Solar: x2.0 QLD1-Solar: x4.0 QLD1-Wind: x4.0 QLD1-Battery: x4.0 SA1-Rooftop Solar: x2.0 SA1-Solar: x4.0 SA1-Wind: x4.0 SA1-Battery: x4.0 TAS1-Battery: x4.0 TAS1-Wind: x4.0 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x2.0 VIC1-Solar: x4.0 VIC1-Wind: x4.0 VIC1-Battery: x4.0	5121.282	48.404	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 12751.0 Hydro: 12808.0 Rooftop Solar: 31550.39 Solar: 61908.0 Wind: 64032.0	123598.332	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
5_5xVreTransitionZeroCoal	NSW1-Black Coal: x0.0 NSW1-Rooftop Solar: x2.0 NSW1-Solar: x5.0 NSW1-Wind: x5.0 NSW1-Battery: x4.0 QLD1-Black Coal: x0.0 QLD1-Rooftop Solar: x2.0 QLD1-Solar: x5.0 QLD1-Wind: x5.0 QLD1-Battery: x4.0 SA1-Rooftop Solar: x2.0 SA1-Solar: x5.0 SA1-Wind: x5.0 SA1-Battery: x4.0 TAS1-Battery: x4.0 TAS1-Wind: x5.0 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x2.0 VIC1-Solar: x5.0 VIC1-Wind: x5.0 VIC1-Battery: x4.0	3622.306	48.404	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 12751.0 Hydro: 12808.0 Rooftop Solar: 31550.39 Solar: 77385.0 Wind: 80040.0	219661.258	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
6.0_6xVreTransitionZeroCoal	NSW1-Black Coal: x0.0 NSW1-Rooftop Solar: x2.0 NSW1-Solar: x6.0 NSW1-Wind: x6.0 NSW1-Battery: x4.0 QLD1-Black Coal: x0.0 QLD1-Rooftop Solar: x2.0 QLD1-Solar: x6.0 QLD1-Wind: x6.0 QLD1-Battery: x4.0 SA1-Rooftop Solar: x2.0 SA1-Solar: x6.0 SA1-Wind: x6.0 SA1-Battery: x4.0 TAS1-Battery: x4.0 TAS1-Wind: x6.0 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x2.0 VIC1-Solar: x6.0 VIC1-Wind: x6.0 VIC1-Battery: x4.0	2856.965	48.404	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 12751.0 Hydro: 12808.0 Rooftop Solar: 31550.39 Solar: 92862.0 Wind: 96048.0	284040.810	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
6.1_6xVreTransitionZeroCoal	NSW1-Black Coal: x0.0 NSW1-Gas: x0.9 NSW1-Rooftop Solar: x2.0 NSW1-Solar: x6.0 NSW1-Wind: x6.0 NSW1-Battery: x4.0 QLD1-Black Coal: x0.0 QLD1-Gas: x0.9 QLD1-Rooftop Solar: x2.0 QLD1-Solar: x6.0 QLD1-Wind: x6.0 QLD1-Battery: x4.0 SA1-Gas: x0.9 SA1-Rooftop Solar: x2.0 SA1-Solar: x6.0 SA1-Wind: x6.0 SA1-Battery: x4.0 TAS1-Gas: x0.9 TAS1-Wind: x6.0 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.0 VIC1-Gas: x0.9 VIC1-Rooftop Solar: x2.0 VIC1-Solar: x6.0 VIC1-Wind: x6.0 VIC1-Battery: x4.0	2697.012	48.404	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 11475.9 Hydro: 12808.0 Rooftop Solar: 31550.39 Solar: 92862.0 Wind: 96048.0	284013.284	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
6.2_6xVreTransitionZeroCoal	NSW1-Black Coal: x0.0 NSW1-Gas: x0.75 NSW1-Rooftop Solar: x2.0 NSW1-Solar: x6.0 NSW1-Wind: x6.0 NSW1-Battery: x4.0 QLD1-Black Coal: x0.0 QLD1-Gas: x0.75 QLD1-Rooftop Solar: x2.0 QLD1-Solar: x6.0 QLD1-Wind: x6.0 QLD1-Battery: x4.0 SA1-Gas: x0.75 SA1-Rooftop Solar: x2.0 SA1-Solar: x6.0 SA1-Wind: x6.0 SA1-Battery: x4.0 TAS1-Gas: x0.75 TAS1-Wind: x6.0 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.0 VIC1-Gas: x0.75 VIC1-Rooftop Solar: x2.0 VIC1-Solar: x6.0 VIC1-Wind: x6.0 VIC1-Battery: x4.0	2441.481	48.404	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 9563.25 Hydro: 12808.0 Rooftop Solar: 31550.39 Solar: 92862.0 Wind: 96048.0	283898.415	NSW1: 0.034 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 11.846
6.3_6xVreTransitionZeroCoal	NSW1-Black Coal: x0.0 NSW1-Gas: x0.5 NSW1-Rooftop Solar: x2.0 NSW1-Solar: x6.0 NSW1-Wind: x6.0 NSW1-Battery: x4.0 QLD1-Black Coal: x0.0 QLD1-Gas: x0.5 QLD1-Rooftop Solar: x2.0 QLD1-Solar: x6.0 QLD1-Wind: x6.0 QLD1-Battery: x4.0 SA1-Gas: x0.5 SA1-Rooftop Solar: x2.0 SA1-Solar: x6.0 SA1-Wind: x6.0 SA1-Battery: x4.0 TAS1-Gas: x0.5 TAS1-Wind: x6.0 TAS1-Battery: x4.0 VIC1-Brown Coal: x0.0 VIC1-Gas: x0.5 VIC1-Rooftop Solar: x2.0 VIC1-Solar: x6.0 VIC1-Wind: x6.0 VIC1-Battery: x4.0	1934.540	48.404	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 6375.5 Hydro: 12808.0 Rooftop Solar: 31550.39 Solar: 92862.0 Wind: 96048.0	283820.702	NSW1: 37.888 QLD1: 40.139 SA1: 3.36 TAS1: 0.0 VIC1: 80.065
7.0_6xVre4xRoofZeroCoal	NSW1-Black Coal: x0.0 NSW1-Gas: x0.5 NSW1-Rooftop Solar: x4.0 NSW1-Solar: x6.0 NSW1-Wind: x6.0 NSW1-Battery: x4.0 QLD1-Black Coal: x0.0 QLD1-Gas: x0.5 QLD1-Rooftop Solar: x4.0 QLD1-Solar: x6.0 QLD1-Wind: x6.0 QLD1-Battery: x4.0 SA1-Gas: x0.5 SA1-Rooftop Solar: x4.0 SA1-Solar: x6.0 SA1-Wind: x6.0 SA1-Battery: x4.0 TAS1-Gas: x0.5 TAS1-Wind: x6.0 TAS1-Battery: x4.0 VIC1-Gas: x0.5 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x4.0 VIC1-Solar: x6.0 VIC1-Wind: x6.0 VIC1-Battery: x4.0	1899.059	48.404	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 6375.5 Hydro: 12808.0 Rooftop Solar: 62848.33 Solar: 92862.0 Wind: 96048.0	335129.580	NSW1: 34.194 QLD1: 39.98 SA1: 2.835 TAS1: 0.0 VIC1: 77.032
8.0_6xVre&BatteryZeroCoal	NSW1-Black Coal: x0.0 NSW1-Gas: x1.0 NSW1-Rooftop Solar: x6.0 NSW1-Solar: x6.0 NSW1-Wind: x6.0 NSW1-Battery: x6.0 QLD1-Black Coal: x0.0 QLD1-Gas: x1.0 QLD1-Rooftop Solar: x6.0 QLD1-Solar: x6.0 QLD1-Wind: x6.0 QLD1-Battery: x6.0 SA1-Gas: x1.0 SA1-Rooftop Solar: x6.0 SA1-Solar: x6.0 SA1-Wind: x6.0 SA1-Battery: x6.0 TAS1-Gas: x1.0 TAS1-Wind: x6.0 TAS1-Battery: x6.0 VIC1-Gas: x1.0 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x6.0 VIC1-Solar: x6.0 VIC1-Wind: x6.0 VIC1-Battery: x6.0	2110.003	72.606	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 12751.0 Hydro: 12808.0 Rooftop Solar: 94146.27 Solar: 92862.0 Wind: 96048.0	327100.374	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
8.1_6xVre&BatteryZeroCoal	NSW1-Black Coal: x0.0 NSW1-Gas: x0.75 NSW1-Rooftop Solar: x6.0 NSW1-Solar: x6.0 NSW1-Wind: x6.0 NSW1-Battery: x6.0 QLD1-Black Coal: x0.0 QLD1-Gas: x0.75 QLD1-Rooftop Solar: x6.0 QLD1-Solar: x6.0 QLD1-Wind: x6.0 QLD1-Battery: x6.0 SA1-Gas: x0.75 SA1-Rooftop Solar: x6.0 SA1-Solar: x6.0 SA1-Wind: x6.0 SA1-Battery: x6.0 TAS1-Gas: x0.75 TAS1-Wind: x6.0 TAS1-Battery: x6.0 VIC1-Gas: x0.75 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x6.0 VIC1-Solar: x6.0 VIC1-Wind: x6.0 VIC1-Battery: x6.0	1701.605	72.606	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 9563.25 Hydro: 12808.0 Rooftop Solar: 94146.27 Solar: 92862.0 Wind: 96048.0	385892.323	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0
8.2_6xVre&BatteryZeroCoal	NSW1-Black Coal: x0.0 NSW1-Gas: x0.5 NSW1-Rooftop Solar: x6.0 NSW1-Solar: x6.0 NSW1-Wind: x6.0 NSW1-Battery: x6.0 QLD1-Black Coal: x0.0 QLD1-Gas: x0.5 QLD1-Rooftop Solar: x6.0 QLD1-Solar: x6.0 QLD1-Wind: x6.0 QLD1-Battery: x6.0 SA1-Gas: x0.5 SA1-Rooftop Solar: x6.0 SA1-Solar: x6.0 SA1-Wind: x6.0 SA1-Battery: x6.0 TAS1-Gas: x0.5 TAS1-Wind: x6.0 TAS1-Battery: x6.0 VIC1-Gas: x0.5 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x6.0 VIC1-Solar: x6.0 VIC1-Wind: x6.0 VIC1-Battery: x6.0	1279.524	72.606	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 6375.5 Hydro: 12808.0 Rooftop Solar: 94146.27 Solar: 92862.0 Wind: 96048.0	385688.088	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 9.694
8.3_VreCurtailReview	NSW1-Black Coal: x0.0 NSW1-Gas: x0.5 NSW1-Rooftop Solar: x3.0 NSW1-Solar: x3.0 NSW1-Wind: x6.0 NSW1-Battery: x8.0 QLD1-Black Coal: x0.0 QLD1-Gas: x0.5 QLD1-Rooftop Solar: x3.0 QLD1-Solar: x3.0 QLD1-Wind: x6.0 QLD1-Battery: x6.0 SA1-Gas: x0.5 SA1-Rooftop Solar: x2.0 SA1-Solar: x3.0 SA1-Wind: x2.0 SA1-Battery: x5.0 TAS1-Gas: x0.5 TAS1-Rooftop Solar: x3.0 TAS1-Wind: x5.0 TAS1-Battery: x6.0 VIC1-Gas: x0.5 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x4.0 VIC1-Solar: x4.0 VIC1-Wind: x4.0 VIC1-Battery: x8.0	1721.569	84.290	Black Coal: 0.0 Brown Coal: 0.0 Diesel: 1110.0 Gas: 6375.5 Hydro: 12808.0 Rooftop Solar: 49482.92 Solar: 48465.0 Wind: 72321.0	176480.790	NSW1: 0.0 QLD1: 0.0 SA1: 0.0 TAS1: 0.0 VIC1: 0.0

In the 0_2024_baseline scenario above, gas-powered generation - GPG (GWh) is much lower than in reality because I’ve set the marginal cost of gas higher than coal and there are no other influencing constraints such as carbon emissions. PyPSA is solving based on least-cost optimisation of power plant and storage dispatch within network constraints. This is much simpler than in real life, which involves a bid stack and many more constraints.

It is worth adding a comparison column to compare more realistic gas generation from Open Electricity to the high VRE scenario that I’m interested in: 8.3_VreCurtailReview . This scenario has zero coal generation and an 84% reduction in gas generation.

The 8.3_VreCurtailReview is one of the more plausible scenarios found by scaling up VRE and storage from the baseline model using a csv import file with different scale-up factors. Please note that, in this scenario, battery storage capacity is 84GW vs the 56 GW forecasted for 2049-50 in the AEMO ISP Step Change scenario for installed capacity.

Show the Python Code

# 1) Filter scenario
df_filt = results_tbl[results_tbl["Scenario"] == "8.3_VreCurtailReview"].copy()

# 2) Compute the % change vs 10,919 GWh. Source: Open Electricity Data (this value will change regularly because it is a rolling 12-month period)
base = 10919
df_filt["Gas Δ% vs 10919GWh"] = (df_filt["GPG (GWh)"] - base) / base * 100

# 3) Insert the new column right after “GPG (GWh)”
col_idx = df_filt.columns.get_loc("GPG (GWh)")
df_filt.insert(col_idx+1, "Gas Δ% vs 10919GWh", df_filt.pop("Gas Δ% vs 10919GWh"))

# Inspect
html_subset_table = df_filt.iloc[:, :5].to_html(escape=False, index=False)
html_subset_table = html_subset_table.replace('<table', 
    '<table style="width: 100%; table-layout: fixed; font-size: 11px; font-family: monospace;"')
html_subset_table = html_subset_table.replace('<td', 
    '<td style="min-width: 250px; word-wrap: break-word; padding: 4px;"')

html_subset_table

Scenario	Objective	GPG (GWh)	Gas Δ% vs 10919GWh	Battery Capacity (GW)
8.3_VreCurtailReview	NSW1-Black Coal: x0.0 NSW1-Gas: x0.5 NSW1-Rooftop Solar: x3.0 NSW1-Solar: x3.0 NSW1-Wind: x6.0 NSW1-Battery: x8.0 QLD1-Black Coal: x0.0 QLD1-Gas: x0.5 QLD1-Rooftop Solar: x3.0 QLD1-Solar: x3.0 QLD1-Wind: x6.0 QLD1-Battery: x6.0 SA1-Gas: x0.5 SA1-Rooftop Solar: x2.0 SA1-Solar: x3.0 SA1-Wind: x2.0 SA1-Battery: x5.0 TAS1-Gas: x0.5 TAS1-Rooftop Solar: x3.0 TAS1-Wind: x5.0 TAS1-Battery: x6.0 VIC1-Gas: x0.5 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x4.0 VIC1-Solar: x4.0 VIC1-Wind: x4.0 VIC1-Battery: x8.0	1721.569	-84.233272	84.29

`8.3_VreCurtailReview` scenario

Load the 8.3_VreCurtailReview network.

Show the Python Code

# Load the network for the specific scenario
n = pypsa.Network("results/scenarios/8.3_VreCurtailReview.nc")

# Metadata to provide as arguments to the plot_dispatch() function.
scenario = "8.3_VreCurtailReview"
objective_text = df_results.loc[df_results["Scenario"] == scenario, "Objective"].values[0]

INFO:pypsa.io:Imported network 8.3_VreCurtailReview.nc has buses, carriers, generators, lines, links, loads, storage_units

Show the Python Code

plot_dispatch(n, time="2024-05-18", days=10, regions=["VIC1"], show_imports=True, scenario_name=scenario, scenario_objective=objective_text, interactive=False)

In the above snapshot there is minimal gas dispatch due to high levels of renewable energy. The trade-off is curtailment, the middle of the day on the 19th of May is a good example, where almost all solar is curtailed.

Show the Python Code

plot_dispatch(n, time="2024-09-20", days=10, regions=["SA1"], show_imports=True, scenario_name=scenario, scenario_objective=objective_text, interactive=False)

As anticipated, curtailment increases dramatically during the shoulder season, as shown above for South Australia.

Show the Python Code

plot_generator_capacity_by_carrier(n)

In this scenario, the wind generation capacity closely matches the 2049-50 ISP Step Change forecast, while solar capacity is lower as a result of my efforts to minimize curtailment. Coal has been reduced to zero and there is only a small amount of Gas. There are no changes to diesel.

Assessing peak demand alongside installed generation capacity for each region in the 8.3_VreCurtailReview network clearly illustrates the substantial overbuilding of renewable energy relative to demand.

Also, consider the relationship between gas generation capacity and demand in SA1, and similarly to an even larger extent, between hydro capacity and demand in TAS1. As a result, these regions tend to experience unserved energy later than other regions.

Show the Python Code

plot_peak_demand_vs_capacity(n, stack_carrier=True)

Total generation per region also highlights the over-generation required compared to baseline.

Show the Python Code

plot_total_demand_vs_generation(n, stack_carrier=True, relative=False)

Show the Python Code

# Note: bus and region are synonymous in this context.
# when calculating curtailment manually (not using inbuild statistics module) we must consider 30-minute timesteps..
def curtailment_region_carrier(n, carriers=['Rooftop Solar', 'Solar', 'Wind']):
    """
    Calculate curtailment % by bus and carrier (excluding Hydro) and pivot carriers wide.
    Adds one monthly curtailed TWh string column per carrier named
    'monthly_curtailment_<carrier>_twh_str'.

    Fixed so that 'month' in the intermediate DataFrame is a Timestamp,
    which lets us do .to_period('M') later without error.
    """
    
    records = []
    monthly_data = {carrier: [] for carrier in carriers}

    # create a PeriodIndex of each snapshot’s month
    snapshot_month = n.snapshots.to_series().dt.to_period('M')

    for carrier in carriers:
        mask = n.generators.carrier == carrier
        gens = n.generators[mask]

        for bus, gens_in_bus in gens.groupby('bus'):
            idx = gens_in_bus.index.intersection(n.generators_t.p_max_pu.columns)

            # overall curtailment % per bus/carrier
            if len(idx) == 0:
                pct = 0.0
            else:
                # Get time step duration automatically
                if len(n.snapshots) > 1:
                    timestep_hours = (n.snapshots[1] - n.snapshots[0]).total_seconds() / 3600
                else:
                    timestep_hours = 1.0  # Default to 1 hour if single snapshot

                pot = (n.generators_t.p_max_pu[idx]
                    .multiply(n.generators.loc[idx, 'p_nom'], axis=1)
                    ).sum().sum() * timestep_hours
                act = n.generators_t.p[idx].sum().sum() * timestep_hours
                curtailed = pot - act
                pct = 100 * curtailed / pot if pot > 0 else 0.0

            records.append({'bus': bus, 'carrier': carrier, 'curtailment_pct': round(pct, 3)})

            # monthly curtailed MWh per bus/carrier
            if len(idx) > 0:
                pot_m = (n.generators_t.p_max_pu[idx]
                        .multiply(n.generators.loc[idx, 'p_nom'], axis=1))
                act_m = n.generators_t.p[idx]

                # Apply timestep correction when summing (reusing timestep_hours from above)
                pot_mon = pot_m.groupby(snapshot_month).sum().sum(axis=1) * timestep_hours
                act_mon = act_m.groupby(snapshot_month).sum().sum(axis=1) * timestep_hours
                curtailed_mon = pot_mon - act_mon
            else:
                # build a zero‐series indexed by each unique period
                curtailed_mon = pd.Series(0.0, index=snapshot_month.unique())

            # *** store month as a Timestamp, not Period ***
            for period, val in curtailed_mon.items():
                monthly_data[carrier].append({
                    'month': period.to_timestamp(),  # <- convert here
                    'bus': bus,
                    'curtailed_mwh': val
                })

    # build the bus×carrier % pivot
    df = pd.DataFrame(records)
    pivot_df = df.pivot(index='bus', columns='carrier', values='curtailment_pct').fillna(0.0)

    # for each carrier, build its monthly TWh‐string column
    for carrier in carriers:
        mon_df = pd.DataFrame(monthly_data[carrier])
        summed = (mon_df.groupby(['bus', 'month'])['curtailed_mwh']
                     .sum()
                     .reset_index())

        # now month is Timestamp, so this works:
        start = summed['month'].min().to_period('M')
        end   = summed['month'].max().to_period('M')
        months_sorted = pd.period_range(start, end, freq='M').to_timestamp()

        ser = {}
        for bus in pivot_df.index:
            subset = summed[summed['bus'] == bus].set_index('month')['curtailed_mwh'].to_dict()
            arr = np.array([round(subset.get(m, 0.0), 3) for m in months_sorted])
            twh = arr / 1e6 # convert MWh to TWh
            ser[bus] = ' '.join(f'{x:.3f}' for x in twh)

        col = f'monthly_curtailment_{carrier.lower().replace(" ", "_")}_twh_str'
        pivot_df[col] = pivot_df.index.map(ser)

    return pivot_df

Leverage the great_tables python library to highlight curtailment data returned in the curtailment_region_carrier function above.

Show the Python Code

pivot_df = curtailment_region_carrier(n)

# from great_tables import GT, md, system_fonts, nanoplot_options

curtailment_tbl = GT(data=pivot_df \
    .reset_index() \
    .round(2) \
    .rename(columns={'bus': 'Region',
                     'Solar': 'Utility Solar',   
                     'monthly_curtailment_rooftop_solar_twh_str': 'Rooftop Solar Curtailment',
                     'monthly_curtailment_solar_twh_str': 'Utility Solar Curtailment',
                     'monthly_curtailment_wind_twh_str': 'Wind Curtailment'
                     })
    )

# Generate great table for curtailment by region and carrier
# Note: Scale-up objective is hardcoded in the source note.
# objective_text.replace('\n', ', ') 
curtailment_gt = curtailment_tbl.tab_header(
        title="NEM Variable Renewable Energy Curtailment by Region"
        ) \
        .tab_spanner(
        label="Curtailment (%)",
        columns=['Rooftop Solar', 'Utility Solar', 'Wind']
        ) \
        .tab_spanner(
            label="Monthly Curtailment Profiles (TWh)",
            columns=['Rooftop Solar Curtailment',
                     'Utility Solar Curtailment',
                     'Wind Curtailment']
        ) \
        .data_color(
        columns=['Rooftop Solar', 'Utility Solar', 'Wind'],
        palette=[
            "#31a354", "#78c679", "#ffffcc",
            "#fafa8c", "#f4cd1e", "#f8910b"],
        domain=[0, 100]
        ) \
        .cols_width(
            {'Wind': '120px', 'Utility Solar': '120px', 'Rooftop Solar': '120px'}
        ) \
        .cols_align(    
             align='center'
        ) \
        .tab_source_note(
            source_note=md(
                "**Calculation**: Curtailment = (Potential energy - Actual energy) / Potential energy * 100."
                                           
            )
        ) \
        .tab_source_note(
            source_note=md(
                "**Note**: there are currently no listed utility scale solar in the TAS1 [Open Electricity Facilities](https://explore.openelectricity.org.au/facilities/tas1/?tech=solar_utility&status=operating,committed)"
                                           
            )
        ) \
        .tab_source_note(
            source_note=md("**Scenario**: 8.3_VreCurtailReview (Zero coal & 84% reduction in GPG).<br><br> **Scale-up objective from 2024 baseline**: *NSW1-Black Coal: x0.0, NSW1-Gas: x0.5, NSW1-Rooftop Solar: x3.0, NSW1-Solar: x3.0, NSW1-Wind: x6.0, NSW1-Battery: x8.0, QLD1-Black Coal: x0.0, QLD1-Gas: x0.5, QLD1-Rooftop Solar: x3.0, QLD1-Solar: x3.0, QLD1-Wind: x6.0, QLD1-Battery: x6.0, SA1-Gas: x0.5, SA1-Rooftop Solar: x2.0, SA1-Solar: x3.0, SA1-Wind: x2.0, SA1-Battery: x5.0, TAS1-Gas: x0.5, TAS1-Rooftop Solar: x3.0, TAS1-Wind: x5.0, TAS1-Battery: x6.0, VIC1-Gas: x0.5, VIC1-Brown Coal: x0.0, VIC1-Rooftop Solar: x4.0, VIC1-Solar: x4.0, VIC1-Wind: x4.0, VIC1-Battery: x8.0*")
        ) \
        .fmt_nanoplot("Rooftop Solar Curtailment",
                      plot_type="bar",
                      options=nanoplot_options(
                          data_bar_fill_color="#FFE066",
                          y_ref_line_fmt_fn="GWh",
                          )
        ) \
        .fmt_nanoplot("Utility Solar Curtailment",
                      plot_type="bar",
                      options=nanoplot_options(
                          data_bar_fill_color="#FDB324",
                          y_ref_line_fmt_fn="GWh",
                          )
        ) \
        .fmt_nanoplot("Wind Curtailment",
                      plot_type="bar",
                      options=nanoplot_options(
                          data_bar_fill_color="#3BBFE5",
                          y_ref_line_fmt_fn="GWh",
                          )
        ) \
        .tab_options(
            source_notes_font_size='x-small',
            source_notes_padding=3,
            heading_subtitle_font_size='small',
            table_font_names=system_fonts("humanist"),
            data_row_padding='1px',
            heading_background_color='#F0FFFF',
            source_notes_background_color='#F0FFF0',
            column_labels_background_color='#F0FFF0',
            table_background_color='snow',
            data_row_padding_horizontal=3,
            column_labels_padding_horizontal=58
        ) \
            .opt_table_outline()

curtailment_gt

NEM Variable Renewable Energy Curtailment by Region
Region	Curtailment (%)			Monthly Curtailment Profiles (TWh)
Region	Rooftop Solar	Utility Solar	Wind	Rooftop Solar Curtailment	Utility Solar Curtailment	Wind Curtailment
NSW1	70.9	60.52	37.6
QLD1	64.18	58.24	36.0
SA1	47.41	60.53	39.85
TAS1	61.98	0.0	43.62
VIC1	64.26	68.16	44.54
Calculation: Curtailment = (Potential energy - Actual energy) / Potential energy * 100.
Note: there are currently no listed utility scale solar in the TAS1 Open Electricity Facilities
Scenario: 8.3_VreCurtailReview (Zero coal & 84% reduction in GPG). Scale-up objective from 2024 baseline: NSW1-Black Coal: x0.0, NSW1-Gas: x0.5, NSW1-Rooftop Solar: x3.0, NSW1-Solar: x3.0, NSW1-Wind: x6.0, NSW1-Battery: x8.0, QLD1-Black Coal: x0.0, QLD1-Gas: x0.5, QLD1-Rooftop Solar: x3.0, QLD1-Solar: x3.0, QLD1-Wind: x6.0, QLD1-Battery: x6.0, SA1-Gas: x0.5, SA1-Rooftop Solar: x2.0, SA1-Solar: x3.0, SA1-Wind: x2.0, SA1-Battery: x5.0, TAS1-Gas: x0.5, TAS1-Rooftop Solar: x3.0, TAS1-Wind: x5.0, TAS1-Battery: x6.0, VIC1-Gas: x0.5, VIC1-Brown Coal: x0.0, VIC1-Rooftop Solar: x4.0, VIC1-Solar: x4.0, VIC1-Wind: x4.0, VIC1-Battery: x8.0

The model applies relatively harsh curtailment to solar and wind, particularly rooftop solar. This is a result of curtailment already being incorporated into utility-scale solar availability (capacity) factors, which I’ve sourced and then calculated from AEMO’s NEMWEB dispatch data.

Almost every large solar farm in Australia’s southeastern states will be forced to switch off at least one-third of the power they generate by 2027 as delays to critical energy transition projects cause major bottlenecks on the electricity grid. - Ryan Cropp, Australian Financial Review, Jul 10, 2025 [link].

Large-scale solar is fast and cost-effective to build, but it’s also the most exposed to curtailment, especially when system security limits or negative prices hit in the middle of a sunny day. - Chris O’Keefe, Clean Energy Council.

The built-in curtailment would likely, at least partially explain why AEMO uses solar and wind traces, rather than solely relying on historical data, in its ISP modelling.

Comparing rooftop solar with utility solar availability (capacity) factors⁴ in our network will help clarify this:

Show the Python Code

solar_cf = n.generators_t.p_max_pu[['SA1-SOLAR', 'SA1-ROOFTOP-SOLAR',
                         'VIC1-SOLAR', 'VIC1-ROOFTOP-SOLAR',
                         'QLD1-SOLAR', 'QLD1-ROOFTOP-SOLAR',
                         'NSW1-SOLAR', 'NSW1-ROOFTOP-SOLAR']].reset_index()

solar_cf['hour'] = solar_cf['snapshot'].dt.hour
hours = sorted(solar_cf['hour'].unique())
even_hours = [h for h in hours if h % 2 == 0]
even_positions = [hours.index(h) for h in even_hours]

combos = [
    ['SA1-SOLAR', 'SA1-ROOFTOP-SOLAR'],
    ['VIC1-SOLAR', 'VIC1-ROOFTOP-SOLAR'],
    ['QLD1-SOLAR', 'QLD1-ROOFTOP-SOLAR'],
    ['NSW1-SOLAR', 'NSW1-ROOFTOP-SOLAR']
]

fig, axes = plt.subplots(
    nrows=2, ncols=2, figsize=(8.5, 8),
    sharex=True, sharey=True   # share y for common scale
)
axes = axes.flatten()

width = 0.35
x = np.arange(len(hours))

for ax, (sol_col, rt_col) in zip(axes, combos):
    sol_data = [solar_cf[solar_cf['hour']==h][sol_col].dropna().values for h in hours]
    rt_data  = [solar_cf[solar_cf['hour']==h][rt_col].dropna().values for h in hours]

    bp1 = ax.boxplot(
        sol_data,
        positions=x - width/2, widths=width, patch_artist=True,
        boxprops=dict(facecolor='#FDB324', edgecolor='black'),
        medianprops=dict(color='blue'),
        whiskerprops=dict(color='black'),
        capprops=dict(color='black'),
        flierprops=dict(markeredgecolor='black')
    )
    bp2 = ax.boxplot(
        rt_data,
        positions=x + width/2, widths=width, patch_artist=True,
        boxprops=dict(facecolor='#FFE066', edgecolor='black'),
        medianprops=dict(color='green'),
        whiskerprops=dict(color='black'),
        capprops=dict(color='black'),
        flierprops=dict(markeredgecolor='black')
    )

    ax.set_title(f"{sol_col} vs {rt_col}")
    ax.set_xticks(even_positions)
    ax.set_xticklabels(even_hours)
    ax.set_facecolor('#F0FFFF')
    # no individual ylabel here

# shared legend
handles = [bp1["boxes"][0], bp2["boxes"][0]]
labels  = ['Utility-scale Solar', 'Rooftop Solar']
fig.legend(handles, labels, loc='lower center', ncol=2, bbox_to_anchor=(0.5, -0.04))

# one common y-axis label
fig.supxlabel('Hour of Day (2024)')
fig.supylabel('Capacity Factor')
fig.suptitle("Utility Solar vs Rooftop Solar Capacity Factors")
# adjust spacing
fig.subplots_adjust(top=0.90, bottom=0.08, left=0.1, right=0.9, hspace=0.2, wspace=0.1)

plt.show()

There is a clear midday plateau in utility solar output (orange boxplots), likely resulting from negative market prices and/or network constraints such as congestion⁵.

In Victoria, the curtailment plateau is especially pronounced. This is likely due to utility-scale solar farms being situated in parts of the network that experience congestion. Congestion arises when the network’s ability to transmit electricity is constrained, requiring operators to curtail solar generation despite sufficient solar resources being available.

Operating utility solar in Victoria (Open Electricity)

Under near-term operating conditions, connection points in Victoria have on average close to 50% curtailment for solar and 30% curtailment for wind. - AEMO, 2025 Enhanced Locational Information (ELI) Report [link]

Further congestion is noted to impact the curtailment at north-west Victoria, due to local transmission network limitations and existing high level of renewable energy generation. - AEMO, 2025 ELI Report.

In contrast, rooftop solar is far more challenging to curtail in real-world conditions. Although mechanisms exist (such as the capability to remotely disconnect rooftop PV via regulations like South Australia’s Smarter Homes program), these are described as a “tool of last resort” and invoked only in extreme security conditions. The technical ability to enforce widespread, routine curtailment of rooftop systems is currently limited; implementing it more broadly would risk significant social license erosion unless there is a massive uptake of consumer energy resources such as batteries and electric vehicles, allowing households to store their own surplus solar for later use.

A similar dynamic is evident with wind output. The dip in wind is too pronounced to be natural drop-off in wind speed during the middle of the day.

Show the Python Code

wind_cf = n.generators_t.p_max_pu[['SA1-WIND', 'SA1-ROOFTOP-SOLAR',
                         'VIC1-WIND', 'VIC1-ROOFTOP-SOLAR',
                         'QLD1-WIND', 'QLD1-ROOFTOP-SOLAR',
                         'NSW1-WIND', 'NSW1-ROOFTOP-SOLAR']].reset_index()

wind_cf['hour'] = wind_cf['snapshot'].dt.hour
hours = sorted(wind_cf['hour'].unique())
even_hours = [h for h in hours if h % 2 == 0]
even_positions = [hours.index(h) for h in even_hours]

combos = [
    ['SA1-WIND', 'SA1-ROOFTOP-SOLAR'],
    ['VIC1-WIND', 'VIC1-ROOFTOP-SOLAR'],
    ['QLD1-WIND', 'QLD1-ROOFTOP-SOLAR'],
    ['NSW1-WIND', 'NSW1-ROOFTOP-SOLAR']
]

fig, axes = plt.subplots(
    nrows=2, ncols=2, figsize=(8.5, 8),
    sharex=True, sharey=True   # share y for common scale
)
axes = axes.flatten()

width = 0.35
x = np.arange(len(hours))

for ax, (sol_col, rt_col) in zip(axes, combos):
    wind_data = [wind_cf[wind_cf['hour']==h][sol_col].dropna().values for h in hours]
    rt_data  = [wind_cf[wind_cf['hour']==h][rt_col].dropna().values for h in hours]

    bp1 = ax.boxplot(
        wind_data,
        positions=x - width/2, widths=width, patch_artist=True,
        boxprops=dict(facecolor='#3BBFE5', edgecolor='black'),
        medianprops=dict(color='blue'),
        whiskerprops=dict(color='black'),
        capprops=dict(color='black'),
        flierprops=dict(markeredgecolor='black')
    )
    bp2 = ax.boxplot(
        rt_data,
        positions=x + width/2, widths=width, patch_artist=True,
        boxprops=dict(facecolor='#FFE066', edgecolor='black'),
        medianprops=dict(color='green'),
        whiskerprops=dict(color='black'),
        capprops=dict(color='black'),
        flierprops=dict(markeredgecolor='black')
    )

    ax.set_title(f"{sol_col} vs {rt_col}")
    ax.set_xticks(even_positions)
    ax.set_xticklabels(even_hours)
    ax.set_facecolor('#F0FFFF')
    # no individual ylabel here

# shared legend
handles = [bp1["boxes"][0], bp2["boxes"][0]]
labels  = ['Wind', 'Rooftop Solar']
fig.legend(handles, labels, loc='lower center', ncol=2, bbox_to_anchor=(0.5, -0.04))

# one common y-axis label
fig.supxlabel('Hour of Day (2024)')
fig.supylabel('Capacity Factor')
fig.suptitle("Wind vs Rooftop Solar Capacity Factors")

fig.subplots_adjust(top=0.90, bottom=0.08, left=0.1, right=0.9, hspace=0.2, wspace=0.1)
plt.show()

05 Conclusion and Next Steps

This simple model suggests that achieving a high-penetration renewable energy grid may be possible, at least in principle. However, there is an inherent trade-off between deploying large amounts of variable renewable energy and managing increased curtailment. The model focuses solely on power flow analysis and does not consider system security or full system cost implications, beyond assigning marginal costs to generators for least-cost dispatch.

Note

There is no assurance the system would operate in a secure state within the 8.3_VreCurtailReview scenario/network.

Incorporating the necessary constraints to ensure secure system operation would introduce significant complexity and is reserved for future analysis. Additionally, this scenario is highly sensitive to increases in demand - a 5% rise, for instance, leads to unserved energy. This highlights the need for substantial overbuilding to create an adequate reserve buffer, which, in turn, further increases curtailment.

During the development of the 8.3_VreCurtailReview scenario, several scenarios with aggressive coal reduction - such as 2_BalancedAggressiveTransition - led to a doubling of gas-powered generation (GPG). This occurred because GPG was required to compensate for the loss of 75% of baseload capacity previously provided by coal. Further details can be found in the Appendix.

The 6.1_6xVreTransitionZeroCoal scenario achieves a 75% reduction in gas generation⁶ however only results in 48.4 GW of battery capacity - less than the 56 GW of installed capacity projected for 2049–50 in the AEMO ISP Step Change scenario⁷. Consequently, this scenario would likely be more economical from a cost perspective compared to 8.3_VreCurtailReview.

The 8.3_VreCurtailReview scenario results in a 46.7% curtailment rate for wind and utility-scale solar, which is more than double the 20% level projected by the ISP for 2050 (2024 ISP, section 4.3, page 52). In addition, since some curtailment is already embedded in the source data for this scenario, the effective curtailment rate is higher again.

A key trade-off emerged when applying different assumptions under the 8.3_VreCurtailReview scenario. Under an economic optimisation approach, curtailment of utility solar and wind reached 46.7%, but storage was used efficiently. By contrast, renewable maximisation reduced curtailment to ~21%, though at the cost of economically inefficient storage utilisation.

Unserved energy predominantly occurred in May and June, particularly the evening of June 13.

Its is also worth re-visiting that 2024 was a challenging year for wind and hydro production, with the presence of numerous high-pressure systems. These systems typically block winter low-pressure weather patterns, leading to reduced wind speeds and rainfall. This makes 2024 a valuable year for stress-testing future scenarios, especially given CSIRO’s projection that such conditions may become more frequent over time.

Mean sea level atmospheric pressure is increasing over Australia, and there has been an increase in the number of high-pressure systems over southern Australia, which bring dry, clear weather and little rainfall. This increase in atmospheric pressure across southern latitudes is a response to climate change [link].

The below plot highlights the relationship between 5-day moving average (MA) Moomba to Adelaide gas pipeline utilisation, barometric pressure, wind speed and temperature at Port Pirie in South Australia.

Moomba to Adelaide gas pipeline utilisation vs barometric pressure, wind speed and temperature at Port Pirie, South Australia.⁸

Interestingly, 2023 also featured blocking high-pressure systems, however was one of the warmest winters on record.

Many of the cold fronts that usually sweep over Australia in winter were also pushed south by blocking high pressure systems in winter. This caused unusually warm and dry weather over parts of southern Australia and played a big role in starving the Alps of snow [link].

Taking a quick look at gas storage levels at the Iona underground storage facility highlights the sensitivity of temperature to gas utilisation in 2023 compared to 2024, which had closer to average winter temperatures.

Gas storage levels at Iona underground storage, near Port Campbell in southwestern Victoria.⁹

Note

I’ve now added a new 5-min resolution network scenario to the web app. The 8.4.1_ZeroThermal scenario is a battery scale-up scenario with no thermal generation. This scenario only covers the problematic months of May and June 2024 where there were multiple low-wind events. The full year at 5-min resolution network is too large to serve via the web app.

Potential future analysis could incorporate:
- increased transmission (lines/links currently stacked on top of each other as they join the same bus/region)
- Explore Power Transfer Distribution Factors (PTDF) - interconnector limits that move with dispatch (like AEMO’s constraint equations)
- accurate geospatial resolution
- remove built-in curtailment from the utility solar and wind data or apply an alternative measurement approach, such as using synthetic wind and solar traces.
- more detailed costs, including capital costs
- carbon emissions.

06 Appendix

Reserve Requirements

As mentioned in the conclusion, there is no certainty around whether or not the 8.3_VreCurtailReview network could operate in a secure state.

It is possible to determine in very rough terms, how often for example the network would be issued Lack of Reserve (LOR)¹⁰ notices by the system operator. Here we’ll use the per region thresholds published by AEMO.

It is worth noting that these reserve requirements are dynamic and always changing, see example below of a recent LOR notice from AEMO in reference to the South Australian region:

LOR notice example showing dynamic reserve requirements in SA1 region

In the count_reserve_breaches_df funtion below, we can pass reserve thresholds per region, published in the 2024 ISP Inputs and Assumptions workbook.xlsx as a rough guide(sheet: Reserves):

{"NSW1":705, "VIC1":550, "QLD1":710, "SA1":195, "TAS1":140}

Show the Python Code

def count_reserve_breaches_df(n, threshold_mw=850, include_storage=True):
    """
    Count how many hours in 2024 the reserve margin
    (available capacity minus load) at each bus dips below threshold_mw.
    Returns a DataFrame with breach counts and breach ratios.

    Parameters
    ----------
    n : pypsa.Network
        A solved PyPSA network.
    threshold_mw : float or dict
        If float, same threshold (in MW) for every bus.
        If dict, keys are bus names and values are thresholds for those buses;
        any bus not in the dict uses the global default (if provided) or 0.
    include_storage : bool, default True
        Whether to include storage units in available capacity calculation.

    Returns
    -------
    pd.DataFrame
        DataFrame with columns:
          'breach_obs': number of hours below threshold
          'breach_ratio': breach_obs / total_snapshots
        Index is bus name.
    """
    # Generator availability
    gen_avail = n.generators_t.p_max_pu.multiply(n.generators.p_nom, axis=1)
    gen_avail_by_bus = gen_avail.T.groupby(n.generators.bus).sum().T
    
    # Storage availability (if requested and storage exists)
    if include_storage and hasattr(n, 'storage_units') and not n.storage_units.empty:
        # Storage can contribute reserves through:
        # 1. Discharge capacity (limited by p_nom and state_of_charge)
        # 2. Charge reduction (if currently charging)
        
        # Maximum discharge capacity
        if hasattr(n, 'storage_units_t') and hasattr(n.storage_units_t, 'p_max_pu') and not n.storage_units_t.p_max_pu.empty:
            # Use time-varying p_max_pu if available
            storage_p_max = n.storage_units_t.p_max_pu.multiply(n.storage_units.p_nom, axis=1)
        else:
            # Fallback: use nominal capacity for all timesteps
            # Create a DataFrame with p_nom values repeated for each timestep
            timesteps = n.snapshots
            n_timesteps = len(timesteps)
            n_storage = len(n.storage_units)
            
            # Create array with p_nom values repeated for all timesteps
            data = np.tile(n.storage_units.p_nom.values, (n_timesteps, 1))
            
            storage_p_max = pd.DataFrame(
                data=data,
                index=timesteps,
                columns=n.storage_units.index
            )
        
        # Account for State of Charge constraints if available
        if hasattr(n.storage_units_t, 'state_of_charge') and not n.storage_units_t.state_of_charge.empty:
            soc = n.storage_units_t.state_of_charge
            # Calculate maximum energy capacity (MWh)
            max_energy = n.storage_units.p_nom * n.storage_units.max_hours
            # Calculate available energy at each timestep (MWh)
            available_energy = soc.multiply(max_energy, axis=1)
            
            # Limit discharge reserves by available energy 
            # For 30-minute reserves (matching snapshot duration):
            energy_limited_discharge = available_energy / 0.5  # MW for 30 minutes
            # Alternative: For 1-hour reserves, use: available_energy / 1.0
            
            # Save original for diagnostic comparison
            storage_p_max_original = storage_p_max.copy()
            
            # Use the minimum of power capacity and energy-limited capacity
            storage_p_max = np.minimum(storage_p_max, energy_limited_discharge)
            
            # Diagnostic: Check which constraint is binding
            power_binding = (storage_p_max.values == storage_p_max_original.values).sum()
            energy_binding = (storage_p_max.values == energy_limited_discharge.values).sum()
            print(f"Times power constraint was binding: {power_binding}")
            print(f"Times energy constraint was binding: {energy_binding}")
        
        # Current storage dispatch (positive = discharge, negative = charge)
        storage_dispatch = n.storage_units_t.p
        
        # Available upward reserve from storage:
        # - If discharging: remaining discharge capacity
        # - If charging: can reduce charging (turn charging into reserve)
        storage_reserve = storage_p_max.subtract(storage_dispatch.abs(), fill_value=0)
        
        # Ensure non-negative reserves
        storage_reserve = storage_reserve.clip(lower=0)
        
        # Group storage reserves by bus
        storage_avail_by_bus = storage_reserve.T.groupby(n.storage_units.bus).sum().T
        
        # Combine generator and storage availability
        all_gen_buses = set(gen_avail_by_bus.columns) if not gen_avail_by_bus.empty else set()
        all_storage_buses = set(storage_avail_by_bus.columns) if not storage_avail_by_bus.empty else set()
        all_avail_buses = sorted(all_gen_buses | all_storage_buses)
        
        gen_avail_by_bus = gen_avail_by_bus.reindex(columns=all_avail_buses, fill_value=0)
        storage_avail_by_bus = storage_avail_by_bus.reindex(columns=all_avail_buses, fill_value=0)
        storage_avail_by_bus = storage_avail_by_bus.reindex(index=gen_avail_by_bus.index, fill_value=0)
        
        total_avail_by_bus = gen_avail_by_bus.add(storage_avail_by_bus, fill_value=0)
    else:
        total_avail_by_bus = gen_avail_by_bus

    # Load data
    loads = n.loads_t.p_set
    load_by_bus = loads.T.groupby(n.loads.bus).sum().T

    # Align all buses
    all_buses = sorted(set(total_avail_by_bus.columns) | set(load_by_bus.columns))
    total_avail_by_bus = total_avail_by_bus.reindex(columns=all_buses, fill_value=0)
    load_by_bus = load_by_bus.reindex(columns=all_buses, fill_value=0)
    total_avail_by_bus = total_avail_by_bus.reindex(index=load_by_bus.index, fill_value=0)

    # Calculate reserve margin
    reserve = total_avail_by_bus.subtract(load_by_bus, fill_value=0)

    # Apply thresholds
    if isinstance(threshold_mw, dict):
        thresh = pd.Series({bus: threshold_mw.get(bus, 0.0) for bus in all_buses})
    else:
        thresh = pd.Series(threshold_mw, index=all_buses)

    # Count breaches
    breaches = (reserve.lt(thresh, axis=1)).sum()
    total_snapshots = reserve.shape[0]

    df = pd.DataFrame({
        'breach_obs': breaches,
        'breach_ratio': breaches / total_snapshots
    })

    return df


# Usage examples:

# For per-region thresholds
# Source: 2024 ISP Inputs and Assumptions workbook.xlsx (sheet: Reserves)
per_region = {"NSW1": 705, "VIC1": 550, "QLD1": 710, "SA1": 195, "TAS1": 140}

# Compare with generators only
results_gen_only = count_reserve_breaches_df(n, threshold_mw=per_region, include_storage=False)

# Include storage (default)
results_with_storage = count_reserve_breaches_df(n, threshold_mw=per_region)

# Print comparison
print("Reserve breach analysis:")
print("\nGenerators only:")
print(results_gen_only)
print("\nWith storage included:")
print(results_with_storage)

Times power constraint was binding: 70499
Times energy constraint was binding: 17341
Reserve breach analysis:

Generators only:
      breach_obs  breach_ratio
bus                           
NSW1        6177      0.351605
QLD1        4308      0.245219
SA1         5258      0.299294
TAS1        3583      0.203950
VIC1        4039      0.229907

With storage included:
      breach_obs  breach_ratio
bus                           
NSW1         818      0.046562
QLD1         616      0.035064
SA1          893      0.050831
TAS1         838      0.047700
VIC1         718      0.040870

New South Wales is attracting the most breaches at just over 37% of observations breached, however when storage is factored in, the number of breaches are dramatically reduced.

We can take a quick look at the gas generation coverage.

Show the Python Code

gas = (
    n.generators
     .query("carrier == 'Gas'")
     .filter(['bus','p_nom'], axis=1)
     .reset_index()
     .rename(columns={'index':'name'})
     )

# add per-region lookup
gas['reserve'] = gas['bus'].map(per_region)

# Reorder columns
gas = gas[['bus', 'p_nom', 'reserve']]

gas

	bus	p_nom	reserve
0	NSW1	1639.0	705
1	QLD1	1601.0	710
2	SA1	1704.0	195
3	TAS1	186.0	140
4	VIC1	1245.5	550

There is enough reserve of gas capacity, however the model will operate under this on many occasions (as highlighted above), unless a custom constraint is added, which will increase solver run times.

VRE Curtailment Stats

Show the Python Code

def calculate_renewable_curtailment_stats(n, technologies=['Wind', 'Solar', 'Rooftop Solar', 'Hydro']):
    """Calculate renewable curtailment statistics"""
    
    # Get curtailment by technology
    curtailment = n.statistics.curtailment().loc[(slice(None), technologies)] / 1e3
    
    # Get actual generation by technology  
    generation = n.statistics.supply(comps=["Generator"]).loc[(slice(None), technologies)] / 1e3
    
    # Calculate results
    results = {}
    
    for tech in technologies:
        tech_curtailment = curtailment.loc[(slice(None), tech)].sum()
        tech_generation = generation.loc[(slice(None), tech)].sum()
        tech_potential = tech_curtailment + tech_generation
        
        results[tech] = {
            'Curtailed (GWh)': tech_curtailment,
            'Generated (GWh)': tech_generation,
            'Potential (GWh)': tech_potential,
            'Curtailment Rate (%)': (tech_curtailment / tech_potential * 100) if tech_potential > 0 else 0
        }
    
    # Overall renewable stats
    total_curtailment = curtailment.sum()
    total_generation = generation.sum()
    total_potential = total_curtailment + total_generation
    
    results['Total'] = {
        'Curtailed (GWh)': total_curtailment,
        'Generated (GWh)': total_generation,
        'Potential (GWh)': total_potential,
        'Curtailment Rate (%)': (total_curtailment / total_potential * 100) if total_potential > 0 else 0
    }
    
    return pd.DataFrame(results).T

# Usage
curtailment_stats = calculate_renewable_curtailment_stats(n, technologies=['Wind', 'Solar'])
curtailment_stats.round(2)

	Curtailed (GWh)	Generated (GWh)	Potential (GWh)	Curtailment Rate (%)
Wind	70718.29	107241.54	177959.83	39.74
Solar	53126.30	34246.55	87372.85	60.80
Total	123844.59	141488.09	265332.68	46.68

Storage v Curtailment Trade-off

The two plots below highlight the different storage behaviour under economic optimisation of storage versus maximising utility solar and wind by reducing curtailment. Curtailment was reduced (from 47.7% to ~21%) by forcing more conservative discharge decisions (by setting a very small negative marginal cost and cost of storage):

Under zero marginal and storage costs:

Storage utilisation statistics under economic optimisation¹¹

n.storage_units["marginal_cost"] = -0.01
n.storage_units["marginal_cost_storage"] = -0.01

After optimisation…

Under -0.01 marginal and storage costs:

Storage utilisation statistics under renewable maximisation

The Perfect Foresight Effect:
Your model knows exactly when high-value discharge opportunities will occur, so it might:
- Skip charging during a renewable surplus if it knows a better opportunity is coming
- Let renewables curtail now to save storage capacity for a more valuable future discharge
- Optimise across the entire time horizon rather than being “greedy” about capturing every bit of renewable energy

Real-World Implications:
This suggests that:
- Sub-optimal storage dispatch (being more conservative) can actually reduce curtailment
- Perfect optimisation might increase curtailment because it’s more selective about when to use storage
- There’s a trade-off between storage utilisation efficiency and renewable energy capture

In practice, system operators often prefer strategies that minimise curtailment (for policy/environmental reasons) even if they’re not perfectly economically optimal.

Claude (Sonnet 4). (2025, August 19). Storage utilisation analysis and renewable curtailment insights. Interactive analysis session. Anthropic.

No Marinus Link

Marinus Link Stage 1, which has recently reached Final Investment Decision (FID), is included in all scenarios, including the baseline. If this link is removed, the VIC1 region experiences 14.3 GWh of unserved energy. Examination of the challenging winter snapshot below highlights that a portion of this shortfall occurs during that period.

Whether investing in a very costly interconnector to address these rare gaps represents a justifiable business case remains open to debate according to some industry pundits.

Unserved energy example in VIC1 with no Marinus Link

Note

On the evening of June 13 at 18:30, wind generation across the entire NEM was 2.084 GW in the 8.3_VreCurtailReview scenario. This scenario assumes between 2–6 times the baseline wind capacity and 5–8 times the baseline battery storage. Despite this, storage resources were fully depleted, leading to unserved energy in the early hours of June 14.

NEM-wide output at 18:30, 13 Jun 2024 - Wind CF: NSW1 0.33%, VIC1 0.62%, SA1 12.5%, TAS1 1.97%

GPG in 2_BalancedAggressiveTransition

The 2_BalancedAggressiveTransition scenario has a 75% reduction in coal and only a mild increase in renewables and storage (see table below).

This indicates that without rapidly expanding renewables and storage, phasing out coal will strain gas generators and related infrastructure, as this scenario sees gas-powered generation rise by 107% compared to Open Electricity consumption over a recent 12 month period.

Show the Python Code

# 1) Filter scenario
df_bat = results_tbl[results_tbl["Scenario"] == "2_BalancedAggressiveTransition"].copy()

# 2) Compute the % change vs 10,919 GWh. Source: Open Electricity Data (note this will change regularly because it is a rolling 12-month period)
base = 10919
df_bat["Gas Δ% vs 10919GWh"] = (df_bat["GPG (GWh)"] - base) / base * 100

# 3) Insert the new column right after “GPG (GWh)”
col_idx = df_bat.columns.get_loc("GPG (GWh)")
df_bat.insert(col_idx+1, "Gas Δ% vs 10919GWh", df_bat.pop("Gas Δ% vs 10919GWh"))

# Inspect
html_bat_table = df_bat.iloc[:, :5].to_html(escape=False, index=False)
html_bat_table = html_bat_table.replace('<table', 
    '<table style="width: 100%; table-layout: fixed; font-size: 11px; font-family: monospace;"')
html_bat_table = html_bat_table.replace('<td', 
    '<td style="min-width: 250px; word-wrap: break-word; padding: 4px;"')

html_bat_table

Scenario	Objective	GPG (GWh)	Gas Δ% vs 10919GWh	Battery Capacity (GW)
2_BalancedAggressiveTransition	NSW1-Black Coal: x0.25 NSW1-Rooftop Solar: x1.25 NSW1-Solar: x1.25 NSW1-Wind: x1.3 NSW1-Battery: x2.0 QLD1-Black Coal: x0.25 QLD1-Rooftop Solar: x1.25 QLD1-Solar: x1.25 QLD1-Wind: x1.3 QLD1-Battery: x2.0 SA1-Rooftop Solar: x1.25 SA1-Solar: x1.25 SA1-Wind: x1.3 SA1-Battery: x2.0 TAS1-Wind: x1.3 TAS1-Battery: x2.0 VIC1-Brown Coal: x0.25 VIC1-Rooftop Solar: x1.25 VIC1-Solar: x1.25 VIC1-Wind: x1.3 VIC1-Battery: x2.0	21090.899	93.157789	24.202

Load the 2_BalancedAggressiveTransition network. This scenario has a 75% reduction in coal and a limited increase in renewables and storage.

Show the Python Code

# Load the network for the specific scenario
n = pypsa.Network("results/scenarios/2_BalancedAggressiveTransition.nc")

# Metadata to provide as arguments to the plot_dispatch() function.
scenario = "2_BalancedAggressiveTransition"
objective_text = df_results.loc[df_results["Scenario"] == scenario, "Objective"].values[0]

INFO:pypsa.io:Imported network 2_BalancedAggressiveTransition.nc has buses, carriers, generators, lines, links, loads, storage_units

Review a previous dispatch period over 10 days in Victoria.

Show the Python Code

plot_dispatch(n, time="2024-05-18", days=10, regions=["VIC1"], show_imports=True, scenario_name=scenario, scenario_objective=objective_text)

GPG in 4_4xVreTransitionZeroCoal

The 4_4xVreTransitionZeroCoal scenario has a zero coal and a 4x increase in wind and utility solar and storage plus a 2x increase in rooftop solar (see table below).

This reduces gas-powered generation by 53% compared to Open Electricity consumption over a recent 12 month period.

Show the Python Code

# 1) Filter scenario
df_4x = results_tbl[results_tbl["Scenario"] == "4_4xVreTransitionZeroCoal"].copy()

# 2) Compute the % change vs 10,919 GWh. Source: Open Electricity Data (note this will change regularly because it is a rolling 12-month period)
base = 10919
df_4x["Gas Δ% vs 10919GWh"] = (df_4x["GPG (GWh)"] - base) / base * 100

# 3) Insert the new column right after “GPG (GWh)”
col_idx = df_4x.columns.get_loc("GPG (GWh)")
df_4x.insert(col_idx+1, "Gas Δ% vs 10919GWh", df_4x.pop("Gas Δ% vs 10919GWh"))

# Inspect
html_4x_table = df_4x.iloc[:, :5].to_html(escape=False, index=False)
html_4x_table = html_4x_table.replace('<table', 
    '<table style="width: 100%; table-layout: fixed; font-size: 11px; font-family: monospace;"')
html_4x_table = html_4x_table.replace('<td', 
    '<td style="min-width: 250px; word-wrap: break-word; padding: 4px;"')

html_4x_table

Scenario	Objective	GPG (GWh)	Gas Δ% vs 10919GWh	Battery Capacity (GW)
4_4xVreTransitionZeroCoal	NSW1-Black Coal: x0.0 NSW1-Rooftop Solar: x2.0 NSW1-Solar: x4.0 NSW1-Wind: x4.0 NSW1-Battery: x4.0 QLD1-Black Coal: x0.0 QLD1-Rooftop Solar: x2.0 QLD1-Solar: x4.0 QLD1-Wind: x4.0 QLD1-Battery: x4.0 SA1-Rooftop Solar: x2.0 SA1-Solar: x4.0 SA1-Wind: x4.0 SA1-Battery: x4.0 TAS1-Battery: x4.0 TAS1-Wind: x4.0 VIC1-Brown Coal: x0.0 VIC1-Rooftop Solar: x2.0 VIC1-Solar: x4.0 VIC1-Wind: x4.0 VIC1-Battery: x4.0	5121.282	-53.097518	48.404

Load the 4_4xVreTransitionZeroCoal network. This scenario has a 100% reduction in coal and a large increase in renewables and storage.

Show the Python Code

# Load the network for the specific scenario
n = pypsa.Network("results/scenarios/4_4xVreTransitionZeroCoal.nc")

# Metadata to provide as arguments to the plot_dispatch() function.
scenario = "4_4xVreTransitionZeroCoal"
objective_text = df_results.loc[df_results["Scenario"] == scenario, "Objective"].values[0]

INFO:pypsa.io:Imported network 4_4xVreTransitionZeroCoal.nc has buses, carriers, generators, lines, links, loads, storage_units

Review the same dispatch period over 10 days in Victoria. Note the emergence of curtailment.

Show the Python Code

plot_dispatch(n, time="2024-05-18", days=10, regions=["VIC1"], show_imports=True, scenario_name=scenario, scenario_objective=objective_text)

Generate Scenarios Function

Below is the function for generating scenarios from a CSV file.

Ensure any amendments are updated in the the input assumptions file and stored in data/inputs/tech_assumptions.csv before execution.

Show the Python Code

# Script to generate a simplified sensitivity analysis template
# that can be plugged into a PyPSA workflow.

###--- GENERATE SCENARIOS ---####
import os

# Intergrate csv input assumptions for tech - this better documents scenarios
def dict_to_multiline_string(d):
    return "\n".join(f"{k}: {v}" for k, v in d.items()) if d else ""

def generate_scenarios(
        network,
        tech_assumptions_path,
        export_dir="results/scenarios"
    ):
    """
    Apply scaling from a technology-assumptions CSV, solve each scenario, export
    the solved network, and return a summary DataFrame.

    Extra metrics collected:
    • **Unserved Energy (GWh)** - total unserved energy across the system.
    • **Generator Capacity (GW)** - total generator capacity by carrier.
      (built-in `n.generators.p_nom`)

    • **GPG (GWh)** - total dispatched energy from carrier 'Gas'.
    • **Wind and Solar Curtailment (GWh)** - renewable curtailment across the whole
      system (built-in `n.statistics.curtailment`).
    • **Battery Capacity (GW)** - total battery storage capacity in the system.

    Parameters
    ----------
    network : pypsa.Network
        The base PyPSA network object to scale up/down.
    tech_assumptions_path : str
        Path to CSV file with scaling assumptions for generators/storage.
    export_dir : str, default "results/scenarios"
        Directory path to export .nc and CSV results.
    

    Returns
    -------
    pandas.DataFrame
        Scenario metrics including battery capacity (GW), generator capacity (GW), curtailment (GWh), unserved energy (GWh) and GPG (GWh).
    """

    
    os.makedirs(export_dir, exist_ok=True)
    results = []

    tech_df = pd.read_csv(tech_assumptions_path)
    scenario_names = tech_df["scenario"].unique()
    timer_all = perf_counter()

    for scenario in scenario_names:
        print(f"\u2192 Scenario: {scenario}")

        # ── copy network ───────────────────────────────────────────────
        n = network.copy(snapshots=network.snapshots)

        # ── apply scaling from CSV ────────────────────────────────────
        df_s = tech_df[tech_df["scenario"] == scenario]
        for _, row in df_s.iterrows():
            component, bus, carrier, scale = (
                row["component"], row["bus"], row["carrier"], row["scale_factor"]
            )
            if component == "generator":
                mask = (n.generators.bus == bus) & (n.generators.carrier == carrier)
                n.generators.loc[mask, "p_nom"] *= scale
                if "p_nom_max" in n.generators.columns:
                    n.generators.loc[mask, "p_nom_max"] *= scale
                n.generators.loc[mask, "p_nom_extendable"] = False
            elif component == "storage_unit":
                mask = (n.storage_units.bus == bus) & (n.storage_units.carrier == carrier)
                n.storage_units.loc[mask, "p_nom"] *= scale
                if "p_nom_max" in n.storage_units.columns:
                    n.storage_units.loc[mask, "p_nom_max"] *= scale
                n.storage_units.loc[mask, "p_nom_extendable"] = False

        # ── solve ──────────────────────────────────────────────────────
        # Change to open source solver if required, default is HiGHS
        # n.optimize()
        n.optimize(solver_name="gurobi")

        # ── unserved energy ───────────────────────────────────────────
        if "Unserved Energy" in n.generators.carrier.values:
            ue_cols = n.generators[n.generators.carrier == "Unserved Energy"].index
            ue_df = n.generators_t.p[ue_cols]
            ue_GWh = ue_df.sum().sum() / 1e3
            ue_by_bus = (
                ue_df.sum(axis=0)
                .groupby(n.generators.loc[ue_cols, "bus"]).sum() / 1e3
            ).round(3).to_dict()
        else:
            ue_GWh, ue_by_bus = 0.0, {}

        # ── gas energy ────────────────────────────────────────────────
        gas_idx = n.generators.index[n.generators.carrier == "Gas"]
        gas_GWh = 0.0
        if len(gas_idx):
            # Get appropriate weights
            weights = n.snapshot_weightings.iloc[:, 0] if n.snapshot_weightings.ndim > 1 else n.snapshot_weightings
            # Multiply power by time weights before summing
            gas_GWh = (n.generators_t.p[gas_idx].multiply(weights, axis=0)).sum().sum() / 1e3
        gas_GWh = round(gas_GWh, 3)
        
        # ── curtailment (built‑in helper) ─────────────────────────
        if hasattr(n, "statistics") and hasattr(n.statistics, "curtailment"):
            curtailment_df = n.statistics.curtailment()
            
            existing_carriers = curtailment_df.index.get_level_values(1).unique()
            wanted_carriers = ['Wind', 'Solar', 'Rooftop Solar']
            valid_carriers = [c for c in wanted_carriers if c in existing_carriers]
            if valid_carriers:
                vre_curt_GWh = curtailment_df.loc[(slice(None), valid_carriers)].sum() / 1e3
            else:
                vre_curt_GWh = 0.0
        else:
            vre_curt_GWh = 0.0
        
        # ── Battery Capacity ─────────────────────────
        battery_carrier = "Battery"  # adjust as needed
        battery_mask = n.storage_units.carrier == battery_carrier

        battery_capacity_GW = 0.0
        if battery_mask.any():
            battery_capacity_GW = n.storage_units.loc[battery_mask, "p_nom"].sum() / 1e3  # MW to GW
            battery_capacity_GW = round(battery_capacity_GW, 3)
        
        # ── Generator capacity by carrier ─────────────────────────
        capacity_by_carrier = (
            n.generators[n.generators.carrier != 'Unserved Energy']
            .groupby("carrier")["p_nom"]
            .sum()
            .round(3)
            .to_dict()
        )
        
        # ── collect metrics ───────────────────────────────────────────
        results.append({
            "Scenario": scenario,
            "Objective": "\n".join(
                f"{r['bus']}-{r['carrier']}: x{r['scale_factor']}" for _, r in df_s.iterrows()
            ),
            # "Total System Cost (B$)": round(n.objective / 1e9, 3),
            # "Total New Capacity (GW)": round(n.generators.p_nom_opt.sum() / 1e3, 3),
            "GPG (GWh)": gas_GWh,
            # "Total Curtailment (GWh)": round(curt_GWh, 3),
            "Battery Capacity (GW)": battery_capacity_GW,
            "Generator Capacity (GW)": dict_to_multiline_string(capacity_by_carrier),
            "Wind & Solar Curtailment (GWh)": round(vre_curt_GWh, 3),
            "Unserved Energy (GWh)": round(ue_GWh, 3),
            "Unserved by Region (GWh)": dict_to_multiline_string(ue_by_bus)
            
        })

        # ── export solved network ─────────────────────────────────────
        n.export_to_netcdf(os.path.join(export_dir, f"{scenario}.nc"))

    df_results = pd.DataFrame(results)
    timestamp = datetime.now().strftime("%Y%m%d_%H%M")
    df_results.to_csv(os.path.join(export_dir, f"scenarios_summary_{timestamp}.csv"), index=False)
    
    total_time = round(perf_counter() - timer_all, 1)
    print(f"All scenarios completed in {total_time} seconds")
    
    return df_results

Run the above function and assign results to df_results for further analysis.

Show the Python Code

# Important: initialise baseline network first to avoid scaling up unwanted network (n)
n = pypsa.Network("results/high-level_nem.nc")
df_results = generate_scenarios(n, tech_assumptions_path="data/inputs/tech_assumptions.csv")

PyPSA Dispatch Streamlit App

App link

Tip: Right-click or Cmd/Ctrl+Click to open in a new tab.

Tip

Like this Quarto document, the Streamlit app is better viewed on a large screen.

Footnotes

Western Australia and the Northern Territory operate separate electricity systems and are not part of the NEM.↩︎
The marginal price for hydro is $8.58 according to the 2024 ISP Inputs and Assumptions workbook; however, I wanted to prevent the solver from underutilising hydro. I also initially set a very small negative marginal price for battery storage to try and make the model behave more conservatively about discharge decisions - however this created unintended consequences such as discharge underutilisation. In real life (no perfect forsight) storage would be retained at a higher state of charge (SOC) as insurance.↩︎
The ramp_limit_up & ramp_limit_down values in generators.csv are slightly different because I’m averaging similar units in NSW1 (NSW1-COAL is a single generator representing all coal generators in region NSW1)↩︎
Rooftop solar capacity factors determined by normalising (dividing each observation) against the region’s maximum output - see import_rooftop_solar.py for python code.↩︎
For detailed information on solar farm capacity factors and why considering them in isolation can be misleading, see this excellent Wattclarity article.↩︎
Compared to 10,919 GWh - recent 12-month Open Electricity ‘1Y-Simplified’ view. Source: Open Electricity ↩︎
2024 ISP, section 6.1, page 66↩︎
Plot created using maps_utilisation_animation.py located under the scripts folder↩︎
Plot created using iona_storage_by_year.py located under the scripts folder↩︎
AEMO LoR factsheet ↩︎
Plots created using storage_utilisation.py located under the scripts folder↩︎

Introduction

Research Objectives

01 Import

02 Create Network

02.01 Load Network Data

02.02 Initialise the network

03 Add Network Components

03 Solve the network

04 Visualisations

0_2024_baseline scenario

8.3_VreCurtailReview scenario

05 Conclusion and Next Steps

06 Appendix

Reserve Requirements

VRE Curtailment Stats

Storage v Curtailment Trade-off

No Marinus Link

GPG in 2_BalancedAggressiveTransition

GPG in 4_4xVreTransitionZeroCoal

Generate Scenarios Function

PyPSA Dispatch Streamlit App

Footnotes

`0_2024_baseline` scenario

`8.3_VreCurtailReview` scenario