Bulk code cleanliness review — sector coupling

[ktehranchi]

Bulk code cleanliness review — sector coupling

Companion to #780, covering the sector-coupling code specifically: heat, transport, industry,
natural gas and building stock — everything behind build_sector.smk and postprocess_sector.smk.

This issue tracks the 24 sector findings that carry behavioural risk. Sector findings that are
mechanically safe are in #782.

24 findings — 3 critical, 11 high, 6 medium, 4 low

Why sector should be fixed first

Scope was derived from the rule graph — each file traced to the .smk rules that invoke it and the modules
that import it. The result is an inversion worth knowing about:

Scope	Files	Lines	Findings	Per 1k	Critical	Crit per 1k	Touch numerics
Whole modeling domain	59	16,557	46	2.78	4	0.24	5 (11%)
Power system model	37	19,509	64	3.28	1	0.05	10 (16%)
Sector coupling	25	14,301	34	2.38	5	0.35	12 (35%)

Sector coupling has the fewest findings and the most dangerous ones: 0.35 criticals per thousand lines
against the power model's 0.05, a sevenfold difference, and it is the only scope where numerics-touching
findings are a third of the total. The power model has the most problems and the fewest consequences; sector is
the reverse.

Critical in this scope

workflow/scripts/eulp.py:95 — Duplicate column names inside _heat_group / _space_heat_group double-count electric and gas heating

out.electricity.heating.energy_consumption.kwh and out.natural_gas.heating.energy_consumption.kwh each appear twice in both _heat_group (lines 61/92, 67/87) and _space_heat_group (lines 97/121, 101/116) because the commercial block was pasted under the residential block without de-duplicating. aggregate_sector filters by membership, not uniqueness, so the selection list keeps the duplicate label and df[['a','a']].sum(axis=1) adds the column twice (verified: pd.DataFrame({'a':[1,2],'b':[10,20]})[['a','a','b']].sum(axis=1) == [12,24]). Every ResStock/ComStock state file that carries those two columns therefore reports roughly 2x the true electric and gas space-heating load, and that inflated space_heating column flows straight into the {end_use}_space-heating_s{simpl}.pkl output of rule build_service_demand (workflow/rules/build_electricity.smk:463). Nothing asserts against it, so the error is invisible.

Evidence

_space_heat_group: ClassVar[list[str]] = [
    # residential
    "out.electricity.heating.energy_consumption.kwh",          # line 97
    ...
    "out.natural_gas.heating.energy_consumption.kwh",          # line 101
    ...
    # commercial
    "out.natural_gas.heating.energy_consumption.kwh",          # line 116
    ...
    "out.electricity.heating.energy_consumption.kwh",          # line 121
]

# and the aggregator, lines 229-231:
        def aggregate_sector(df: pd.DataFrame, columns: list[str]) -> pd.Series:
            sector_columns = [x for x in columns if x in df.columns.to_list()]
            return df[sector_columns].sum(axis=1)

Fix: De-duplicate while preserving order at class-definition time, and make the invariant enforceable rather than hoping reviewers spot a repeated string in a 30-line literal:

@staticmethod
def _unique(cols: list[str]) -> list[str]:
    return list(dict.fromkeys(cols))

def _aggregate_data(self, df):
    def aggregate_sector(df, columns):
        sector_columns = self._unique([x for x in columns if x in df.columns])
        return df[sector_columns].sum(axis=1)

Better still, define the residential and commercial column sets as separate constants and build the groups as explicit set unions (_RES_SPACE_HEAT | _COM_SPACE_HEAT), which makes the collision structurally impossible. Add a module-level test asserting len(g) == len(set(g)) for every ClassVar group. Note the same paste also produced the dead entry "out.electricity.cooling.energy_consumption.kwh," (line 144, trailing comma inside the string) which silently never matches.

_{Risk: numerics — needs a validation run}

workflow/scripts/plot_sankey_carbon.py:217 — Carbon Sankey divides CO2 tonnes by the energy constant TBTU_2_MWH while labelling the chart "MT"

Every carbon Sankey link value is scaled by ~1/293,000 and then presented with a mass unit suffix, so the published carbon flow charts are numerically meaningless and there is no error, warning, or visual cue. valueformat=".0f" (line 273) then rounds most links to 0. This is the classic copy-paste-then-forget-to-adapt failure, and the surrounding block confirms the provenance: the comment at plot_sankey_carbon.py:17-18 still says "Energy Services and Rejected Energy links do not follow node color assignment", and the mock entrypoint at :228 still names the wrong rule, "plot_sankey_energy".

Evidence

plot_sankey_carbon.py:216-217:

# convert units
df["value"] = df.value.mul(1 / TBTU_2_MWH)  # MWH -> TBTU

The values are CO2 mass drawn off the co2 bus — plot_sankey_carbon.py:101-111:

links = (
    n.links_t["p2"][links.index]
    .mul(-1)
    .mul(weights, axis=0)

and the figure is labelled as mass, plot_sankey_carbon.py:274 and :296:

                valuesuffix="MT",

...
title_text=f"{state} Carbon Flow in {investment_period} (MT)",

The line is a verbatim copy from plot_sankey_energy.py:591, where the same conversion is correct (that chart is `valuesuffix="TBTU"`). constants.py:39: `TBTU_2_MWH = 1e6 * (1 / MMBTU_MWHthemal)` ≈ 2.93e5.

Fix: Remove the unit conversion from plot_sankey_carbon.py's format_sankey_data (or parameterise it, since that function is otherwise 32/39 lines identical to the energy version — see the duplication finding). Confirm the desired output unit against n.links_t.p2 on the co2 bus (tonnes per the network schema) and scale explicitly, e.g. .mul(1e-6) for tonnes -> MT. Also fix the stale rule name at :228 and the stale COLORS comment at :17-18.

_{Risk: numerics — needs a validation run}

workflow/scripts/plot_validation_sector.py:1 — Entire 758-line module is unreachable; the rule that should run it points at a script that crashes

758 lines — 10% of the domain — are maintained by nobody and read by nothing, while the rule intended to use them is broken. plot_sector_validation is also not requested by sector_figures in workflow/Snakefile:200-245, so the breakage is latent: anyone who wires the validation figures into rule all gets a KeyError on params["result"] rather than plots. Meanwhile the module quietly duplicates seven plot functions, a PlottingData dataclass, SECTOR_MAPPER, FIG_WIDTH/FIG_HEIGHT/EXT, and a byte-identical save_fig from plot_statistics_sector.py, so any fix to the live sector plots has a silent twin waiting to diverge.

Evidence

plot_validation_sector.py:651-663 (its only entrypoint):

if name == "main":
if "snakemake" not in globals():
snakemake = mock_snakemake(
"plot_sector_validation",

But workflow/rules/postprocess_sector.smk:172-194 declares:

rule plot_sector_validation:
params:
plotting=config["plotting"],
eia_api=config["api"]["eia"],
root_dir=RESULTS
+ "{interconnect}/figures/s{simpl}c{clusters}/l{ll}{opts}_{sector}/",
...
script:
"../scripts/plot_statistics_sector.py"

and plot_statistics_sector.py:1525-1526 does, unconditionally:

params = snakemake.params
result = params["result"]

`grep -rn "plot_validation_sector" .` returns zero hits repo-wide. Every other sector plotting rule (postprocess_sector.smk:60, :84, :108) passes `result="emissions"/"production"/"capacity"`; `plot_sector_validation` is the only one that does not.

Fix: Decide the module's fate explicitly. If sector validation is wanted: point rule plot_sector_validation at ../scripts/plot_validation_sector.py, and reconcile FIGURES_SECTOR_VALIDATE = ["emissions_by_sector", ...] (postprocess_sector.smk:25-30) with VALIDATION_PLOTS names in plot_validation_sector.py:570-619, which currently emit emissions_by_sector_validation — the names do not match today either. If it is not wanted, delete plot_validation_sector.py and rule plot_sector_validation together, and drop FIGURES_SECTOR_VALIDATE.

_{Risk: numerics — needs a validation run}

Two further sector criticals are one-line, no-behaviour-change fixes tracked in #782:
plot_statistics_sector.py:200 (three plot functions filter on a module-global month_i instead of their own
month parameter, so a correctly-labelled chart can show the wrong month) and build_co2_storage.py:10.

---

Findings

Error handling (8)

• workflow/scripts/add_sectors.py:87 — co2_carrier defaults to the literal string "carrier", and the resulting KeyError is silently swallowed
  _{high · Error handling · risk: numerics}
• workflow/scripts/add_sectors.py:749 — scale_exising_stock is read in the industrial branch but only bound inside the service-sector branch
  _{high · Error handling · risk: low}
• workflow/scripts/build_stock_data.py:1448 — Water-heater brownfield called with fuel="oil" but dispatch only matches "elec"/"gas"/"lpg"
  _{high · Error handling · risk: low}
• workflow/scripts/plot_statistics_sector.py:153 — 24 blanket except TypeError handlers wrap 5-8 statements each and turn any failure into a blank figure
  _{high · Error handling · risk: low}
• workflow/scripts/build_transportation.py:418 — isinstance(marginal_cost, None) raises TypeError instead of handling the None case it was written for
  _{medium · Error handling · risk: low}
• workflow/scripts/opts/sector.py:342 — Bare except KeyError: continue with the warning commented out silently drops gas trade limits
  _{medium · Error handling · risk: low}
• workflow/scripts/opts/sector.py:163 — Missing co2 policy config logs an error and returns, leaving an uncapped model that still solves
  _{low · Error handling · risk: low}
• workflow/scripts/summary_sector.py:112 — _resample_data silently returns unresampled data; a bare string sits where a raise belongs
  _{low · Error handling · risk: low}

Duplication (6)

• workflow/scripts/eulp.py:95 — Duplicate column names inside _heat_group / _space_heat_group double-count electric and gas heating
  _{CRITICAL · Duplication · risk: numerics}
• workflow/scripts/plot_sankey_carbon.py:217 — Carbon Sankey divides CO2 tonnes by the energy constant TBTU_2_MWH while labelling the chart "MT"
  _{CRITICAL · Duplication · risk: numerics}
• workflow/scripts/build_heat.py:524 — Four near-identical load-splitting copies, all misusing str.rstrip as a suffix strip
  _{high · Duplication · risk: low}
• workflow/scripts/plot_statistics_sector.py:479 — Three capacity plotters ignore period in their loop; every horizon subplot renders identical data
  _{high · Duplication · risk: numerics}
• workflow/scripts/build_emission_tracking.py:113 — _add_ch4_carrier reads the CO2 tech colour, so methane plots in carbon dioxide's colour
  _{low · Duplication · risk: low}
• workflow/scripts/summary_sector.py:331 — Triplicated capacity-filter preamble has diverged: only the brownfield copy special-cases 'trn'
  _{low · Duplication · risk: numerics}

Testability (3)

• workflow/scripts/add_extra_components.py:1600 — costs leaks out of the horizon loop; PTC and CO2 logic silently use the last horizon's cost table
  _{high · Testability · risk: numerics}
• workflow/scripts/build_transportation.py:439 — Zero tests for sector coupling; entry points need a full network, disk paths and a live EIA key
  _{medium · Testability · risk: low}
• workflow/scripts/opts/sector.py:422 — add_ng_import_export_limits performs live EIA HTTP requests from inside extra_functionality
  _{medium · Testability · risk: low}

Nesting (2)

• workflow/scripts/plot_statistics_sector.py:1646 — Per-state monthly figure loop sits outside the state loop, so only the last state gets month plots
  _{high · Nesting · risk: numerics}
• workflow/scripts/add_extra_components.py:1139 — add_co2_storage: 196 lines building two unrelated topologies behind one five-deep if sector: pyramid
  _{medium · Nesting · risk: low}

Dead code (2)

• workflow/scripts/plot_validation_sector.py:1 — Entire 758-line module is unreachable; the rule that should run it points at a script that crashes
  _{CRITICAL · Dead code · risk: numerics}
• workflow/scripts/build_demand.py:617 — ReadEulp._format_data throws away _apply_timeshift, so EULP profiles are never timezone-shifted
  _{high · Dead code · risk: numerics}

Function size (1)

• workflow/scripts/build_stock_data.py:973 — add_service_brownfield is a 481-line function with 8 nested closures (cyclomatic complexity 60)
  _{high · Function size · risk: numerics}

Naming (1)

• workflow/scripts/opts/sector.py:477 — Two nested add_capacity_constraint functions in one module read shift in different units
  _{high · Naming · risk: numerics}

Config coupling (1)

• workflow/scripts/summary_sector.py:847 — Transport historical getter hardcodes 2020 while every sibling takes the investment year
  _{medium · Config coupling · risk: numerics}

---

How this list was produced

A ruff sweep with an expanded rule set (F,C90,ARG,ERA,SIM,RET,PLR,B,PIE,TRY,PTH,DTZ) established an
objective baseline, then twelve reviewers each read one domain of the codebase in full rather than sampling.
Eight of the twelve domains went through an adversarial verification pass instructed to refute: open the file
at the line, confirm the quoted evidence is verbatim, grep the whole repo before accepting any dead-code claim,
and diff both sites before accepting any duplication claim. That pass adjusted 64 findings and confirmed 32.
The remaining four domains self-verified under the same instructions, and the critical findings were re-checked
by hand.

Every item below was confirmed by opening the file at the line. Line numbers are against develop as of
2026-08-24. No code has been changed — this is a survey only.

Each finding carries a risk label:

• safe — provably no behaviour change (mostly deletions)
• low — small blast radius, no numerics
• numerics — touches model results; should not land without a validation run

_{Filed from an automated clean-code review. Items are unchecked survey results, not agreed work — triage
before acting. Ping if you want any line promoted into its own issue.}

[ktehranchi]

@trevorb1 I asked claude to do a code cleanliness review on the entire repo, broke out the sector related items here in case you wanted to look. No need to or rush, especially since your defense is coming up! good luck with your prep for that!!