Capacity expansion with UCBlock¶

This notebook depicts the current status of the conversion between SMS++ and PyPSA using the following tools under development:

• pySMSpp: it aims to provide an abstract python interface for SMS++ models
• pypsa2smspp: it aims to provide a converter from PyPSA to SMS++ models exploiting pySMSpp

In this notebook, the following steps are performed:

1. Create a simple PyPSA model
2. Convert the PyPSA model to SMS++ using pypsa2smspp and optimize the model using default conversion to SMS++ UCBlock component
3. Verify the equivalence of objective function
4. Analyze the PyPSA model created with SMS++ and pypsa2smspp and verify dispatch

This notebook will use the UCBlock solver from SMS++ tools to optimize the model. To do so, pypsa2smspp will convert the PyPSA model to an UCBlock and optimize it using the default settings.

1. Creation of the PyPSA model¶

We create a simple PyPSA model arbitrary number of buses, few generators, storages and loads.

Main assumptions are:

• The network is purely radial in the form bus1 -> bus2 -> bus3 -> ... busN
• A load is added to each bus of the network
• A generator is added to the first node

In [1]:

import pypsa
import pandas as pd
import pypsa2smspp

import pypsa import pandas as pd import pypsa2smspp

The following code creates the desired pypsa model.¶

In [2]:

n = pypsa.Network()
n.set_snapshots(pd.date_range("2024-01-01T00:00", "2024-01-01T23:00", freq="h"))

# Add carriers
n.add("Carrier", "AC")

# Add buses
n_buses = 2
for b in range(n_buses):
    n.add("Bus", f"bus{b}", carrier="AC")

# Add lines in a radial topology using bidirectional links
n_lines = n_buses - 1
for l in range(n_lines):
    n.add(
        "Link",
        f"line{l}",
        bus0=f"bus{l}",
        bus1=f"bus{l+1}",
        length=1,
        capital_cost=1000,
        p_min_pu=-1,
        p_nom_extendable=True,
    )

# Add a load to each bus
n_loads = n_buses
for l in range(n_loads):
    n.add("Load", f"load{l}", bus=f"bus{l}", p_set=pd.Series(100, index=n.snapshots))

# Add a generator to the first bus
n.add(
    "Generator",
    "gen0",
    bus="bus0",
    p_nom_extendable=True,
    capital_cost=1000,
    marginal_cost=1,
)

n = pypsa.Network() n.set_snapshots(pd.date_range("2024-01-01T00:00", "2024-01-01T23:00", freq="h")) # Add carriers n.add("Carrier", "AC") # Add buses n_buses = 2 for b in range(n_buses): n.add("Bus", f"bus{b}", carrier="AC") # Add lines in a radial topology using bidirectional links n_lines = n_buses - 1 for l in range(n_lines): n.add( "Link", f"line{l}", bus0=f"bus{l}", bus1=f"bus{l+1}", length=1, capital_cost=1000, p_min_pu=-1, p_nom_extendable=True, ) # Add a load to each bus n_loads = n_buses for l in range(n_loads): n.add("Load", f"load{l}", bus=f"bus{l}", p_set=pd.Series(100, index=n.snapshots)) # Add a generator to the first bus n.add( "Generator", "gen0", bus="bus0", p_nom_extendable=True, capital_cost=1000, marginal_cost=1, )

In [3]:

n.optimize(solver_name="highs")

n.optimize(solver_name="highs")

/tmp/ipykernel_1813/3144390136.py:1: FutureWarning: The default value of `include_objective_constant` will change from True to False in version 2.0. Set `include_objective_constant` explicitly to suppress this warning. Using False improves LP numerical conditioning by not including the objective constant as a variable.
  n.optimize(solver_name="highs")

INFO:linopy.model: Solve problem using Highs solver

INFO:linopy.io: Writing time: 0.04s

INFO:linopy.constants: Optimization successful: 
Status: ok
Termination condition: optimal
Solution: 50 primals, 146 duals
Objective: 3.05e+05
Solver model: available
Solver message: Optimal

INFO:pypsa.optimization.optimize:The shadow-prices of the constraints Generator-ext-p-lower, Generator-ext-p-upper, Link-ext-p-lower, Link-ext-p-upper were not assigned to the network.

Running HiGHS 1.11.0 (git hash: 364c83a): Copyright (c) 2025 HiGHS under MIT licence terms
LP   linopy-problem-4_t760zo has 146 rows; 50 cols; 242 nonzeros
Coefficient ranges:
  Matrix [1e+00, 1e+00]
  Cost   [1e+00, 1e+03]
  Bound  [0e+00, 0e+00]
  RHS    [1e+02, 1e+02]
Presolving model
72 rows, 2 cols, 72 nonzeros  0s
0 rows, 0 cols, 0 nonzeros  0s
Presolve : Reductions: rows 0(-146); columns 0(-50); elements 0(-242) - Reduced to empty
Solving the original LP from the solution after postsolve
Model name          : linopy-problem-4_t760zo
Model status        : Optimal
Objective value     :  3.0480000000e+05
P-D objective error :  0.0000000000e+00
HiGHS run time      :          0.00

Out[3]:

('ok', 'optimal')

2. Convert the PyPSA model to SMS++ using pypsa2smspp and solve it¶

In [4]:

tran = pypsa2smspp.Transformation()
n_smspp = tran.run(n, verbose=False)

tran = pypsa2smspp.Transformation() n_smspp = tran.run(n, verbose=False)

Executing command:
ucblock_solver temp_test_case.nc -S uc_solverconfig.txt -c /home/docs/checkouts/readthedocs.org/user_builds/pypsa2smspp/conda/latest/lib/python3.13/site-packages/pysmspp/data/configs/UCBlock/ -p /home/docs/checkouts/readthedocs.org/user_builds/pypsa2smspp/checkouts/latest/docs/examples/output/ -O /home/docs/checkouts/readthedocs.org/user_builds/pypsa2smspp/checkouts/latest/docs/examples/output/solution_test_case.nc

Using a default Block configuration
Running HiGHS 1.14.0 (git hash: n/a): Copyright (c) 2026 under MIT licence terms
Solver: HiGHSMILPSolver
LP has 240 rows; 122 cols; 432 nonzeros
Coefficient ranges:
  Matrix  [1e+00, 1e+00]
  Cost    [1e+00, 1e+03]
  Bound   [1e+07, 1e+09]
  RHS     [1e+02, 1e+02]
WARNING: Problem has some excessively large column bounds
WARNING:    Consider scaling the    bounds by 1e-3, or setting the user_bound_scale option to -10
Presolving model
0 rows, 0 cols, 0 nonzeros 0s
0 rows, 0 cols, 0 nonzeros 0s
Presolve reductions: rows 0(-240); columns 0(-122); nonzeros 0(-432) - Reduced to empty
Performed postsolve
Solving the original LP from the solution after postsolve

Model status        : Optimal
Objective value     :  3.0480000000e+05
P-D objective error :  0.0000000000e+00
HiGHS run time      :          0.00
Elapsed time: 9.82821000e-04 s
Status = 10 (Success)
Upper bound = 3.04800000e+05
Lower bound = 3.04800000e+05

----- IntermittentUnitBlock 0 -----
Function value = 2.04800000e+05
Capacity       = 2.00000000e+11
Active power   = [      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02      2.00000000e+02 ]

----- DesignNetworkBlock 0 -----
Function value   = 1.00000000e+05
Design variables = [
  line    0 : 100
]

  --- DCNetworkBlock 0 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 1 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 2 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 3 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 4 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 5 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 6 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 7 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 8 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 9 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 10 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 11 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 12 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 13 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 14 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 15 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 16 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 17 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 18 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 19 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 20 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 21 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 22 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

  --- DCNetworkBlock 23 ---
  Power flow       = [      1.00000000e+02 ]
  Auxiliary variable     = [      1.00000000e+02 ]

Peak CPU Usage: 0.00 %
Peak Memory Usage: 0.92 MB
Total Time: 0.20 seconds

[split] No merged DCNetworkBlock entries found; nothing to do.

Objective value

In [5]:

n_smspp.objective

n_smspp.objective

Out[5]:

304800.0

3. Verify the equivalence of SMS++ and PyPSA results¶

In [6]:

pypsa_obj = n.objective
smspp_obj = n_smspp.objective
print("SMS++ obj         : %.6f" % smspp_obj)
print("PyPSA obj         : %.6f" % pypsa_obj)
print("Error SMS++ - PyPSA [%%]: %.5f" % (100*(smspp_obj - pypsa_obj)/pypsa_obj))

pypsa_obj = n.objective smspp_obj = n_smspp.objective print("SMS++ obj : %.6f" % smspp_obj) print("PyPSA obj : %.6f" % pypsa_obj) print("Error SMS++ - PyPSA [%%]: %.5f" % (100*(smspp_obj - pypsa_obj)/pypsa_obj))

SMS++ obj         : 304800.000000
PyPSA obj         : 304800.000000
Error SMS++ - PyPSA [%]: 0.00000

4. Analyze the PyPSA model created with SMS++ and pypsa2smspp¶

Retrieve PyPSA statistics

In [7]:

n_stat = n.statistics.energy_balance(comps=["Generator"]).droplevel([0, 2])
n_stat

n_stat = n.statistics.energy_balance(comps=["Generator"]).droplevel([0, 2]) n_stat

Out[7]:

carrier
-    4800.0
dtype: float64

Retrieve SMS++ statistics

In [8]:

n_parsed_stat = n_smspp.statistics.energy_balance(comps=["Generator"]).droplevel([0, 2])
n_parsed_stat

n_parsed_stat = n_smspp.statistics.energy_balance(comps=["Generator"]).droplevel([0, 2]) n_parsed_stat

Out[8]:

carrier
-    4800.0
dtype: float64

Calculate difference between the two

In [9]:

error_stat = n_stat - n_parsed_stat

merged_stat = pd.concat([n_stat, n_parsed_stat, error_stat], axis=1)
merged_stat.columns = ["pypsa", "smspp", "difference"]
merged_stat

error_stat = n_stat - n_parsed_stat merged_stat = pd.concat([n_stat, n_parsed_stat, error_stat], axis=1) merged_stat.columns = ["pypsa", "smspp", "difference"] merged_stat

Out[9]:

	pypsa	smspp	difference
carrier
-	4800.0	4800.0	0.0