Resiliency Module

API reference for sdom.resiliency – the operational-resiliency evaluation module. See the user-guide page Resiliency Evaluation for the mathematical background, file-format requirements and worked examples.

Top-Level Convenience

sdom.resiliency.evaluate_resiliency(snapshot_dir, *, inputs_dir, outage_spec, year=2030, scenario_id=1, n_hours=8760, hours=None, min_soc_per_tech=None, slack_penalty=10000.0, curtailment_penalty=0.0, soc_slack_penalty=1000.0, formulation_overrides=None, n_workers=None, solver='highs', solver_options=None, critical_load_MW=None, profile_baseline=False, profile_outages=False)[source]

End-to-end helper: load -> baseline dispatch -> outage evaluation.

Parameters:

  • snapshot_dir (str or pathlib.Path) – Snapshot directory passed to load_designed_system().

  • inputs_dir (str or pathlib.Path) – Previous-stage inputs directory.

  • outage_spec (OutageSpec) – Outage scenario specification.

  • year (int, optional) – Calendar year of the snapshot. Default 2030.

  • scenario_id (int, optional) – Scenario / Run id resolved from the snapshot CSVs. Default 1.

  • n_hours (int, optional) – Baseline-dispatch horizon length. Default 8760.

  • hours (iterable of int, optional) – Anchor hours to evaluate. None (default) evaluates every hour 1..n_hours.

  • min_soc_per_tech (dict, optional) – Operational SOC floor per storage tech (fraction of Cap_E); forwarded to both the baseline builder and the per-hour outage runner. Default None.

  • slack_penalty (float, optional) – Penalty (USD/MWh) applied to slack u[t] in the outage LP. Default 10_000.

  • curtailment_penalty (float, optional) – Penalty applied to curtailed VRE energy (USD/MWh). Default 0.

  • soc_slack_penalty (float, optional) – Penalty (USD/MWh) on the per-tech SOC recovery-target slack in the outage LP. Forwarded to run_resiliency_evaluation() and build_outage_dispatch(). Default 1_000.

  • formulation_overrides (dict, optional) – Component formulation overrides forwarded to load_designed_system().

  • n_workers (int, optional) – Worker pool size for the per-hour evaluation. None (default) resolves to max(1, os.cpu_count() - 1) inside run_resiliency_evaluation().

  • solver (str, optional) – Pyomo solver name. "highs" first tries appsi_highs. Default "highs".

  • solver_options (dict, optional) – Solver options forwarded to both the baseline solve and every per-hour outage solve.

  • critical_load_MW (float, optional) – Constant critical load (MW) used in place of the hourly load series during the outage sub-horizon of each per-hour LP. Forwarded to run_resiliency_evaluation() and through it to build_outage_dispatch(). None (default) preserves the original behaviour. Must be non-negative.

  • profile_baseline (bool, optional) – When True, attach a ModelInitProfiler to the baseline build/solve and print summary tables. Default False (opt-in to avoid runtime/logging overhead in the top-level helper).

  • profile_outages (bool, optional) – When True and n_workers == 1, profile every per-hour outage build. Ignored when n_workers > 1 (a per-worker summary is rarely useful). Default False (opt-in).

Returns:

Per-hour records (sorted by anchor hour) plus run metadata including n_workers_used, n_hours, solver and the outage_spec reference.

Return type:

ResiliencyResults

See also

sdom.resiliency.load_designed_system, sdom.resiliency.build_baseline_dispatch, sdom.resiliency.run_baseline_dispatch, sdom.resiliency.run_resiliency_evaluation

Examples

>>> from sdom.resiliency import OutageSpec, evaluate_resiliency
>>> spec = OutageSpec(
...     duration_hours=4,
...     recovery_hours=4,
...     outaged_assets={"imports": "all"},
... )
>>> results = evaluate_resiliency(
...     "snapshot/",
...     inputs_dir="inputs/",
...     outage_spec=spec,
...     n_hours=24,
...     hours=[1, 5, 10],
...     n_workers=1,
... )

Data Loading

sdom.resiliency.load_designed_system(snapshot_dir, *, inputs_dir, year=2030, scenario_id=1, formulation_overrides=None, attach_cem_data=True)[source]

Load a fixed-capacity designed system from SDOM output snapshots.

Parameters:

  • snapshot_dir (str or pathlib.Path) – Directory containing OutputSummary_* and OutputSelectedVRE_* CSVs from a prior SDOM design run. OutputGeneration_* and OutputStorage_* files are ignored at this stage; files containing Phase1 in their names are also excluded.

  • inputs_dir (str or pathlib.Path) – Directory holding the previous-stage time-series and parameter CSVs (Load_hourly_*, CFSolar_*, CFWind_*, Data_Balancing_units_*, StorageData_*, Import_Cap_*, import_prices_*, Export_Cap_*, export_prices_*, fixed_dem_charges.csv, var_dem_charges.csv).

  • year (int, optional) – Calendar year used in CSV filenames. Default 2030.

  • scenario_id (int, optional) – Scenario / Run id to extract from snapshot CSVs. Default 1. If a snapshot file contains a single unique scenario, that one is used and a warning is emitted when scenario_id differs. If multiple scenarios are present, scenario_id must match one of them; otherwise ValueError is raised listing the available ids.

  • formulation_overrides (dict, optional) – Mapping {component: formulation_name} overlaid on the default formulation map.

  • attach_cem_data (bool, optional) – When True (default), also load the CEM-shaped data dict via load_cem_data() and attach it as DesignedSystem.cem_data. Required by sdom.resiliency.build_baseline_dispatch(). Set to False to skip this step in environments where the previous-stage inputs directory is not a valid CEM inputs folder.

Return type:

DesignedSystem

Raises:

  • FileNotFoundError – When the OutputSummary_* snapshot is missing.

  • ValueError – When scenario_id cannot be resolved against the snapshot.

Examples

>>> ds = load_designed_system(
...     "Data/resiliency_eval/3MW_critical_load_24hrs_outage_24hrs_recovery",
...     inputs_dir="Data/resiliency_eval/inputs_previous_stage/Paper_PGnE/Paper",
...     year=2030,
...     scenario_id=1,
... )
>>> ds.scenario_id
1
class sdom.resiliency.DesignedSystem(storage_caps: dict[str, dict[str, float]] = <factory>, thermal_caps: dict[str, dict[str, float]] = <factory>, solar_caps: dict[str, float] = <factory>, wind_caps: dict[str, float] = <factory>, solar_fom: dict[str, float] = <factory>, wind_fom: dict[str, float] = <factory>, load: ~pandas.core.series.Series | None = None, cf_solar: ~pandas.core.frame.DataFrame | None = None, cf_wind: ~pandas.core.frame.DataFrame | None = None, nuclear: ~pandas.core.series.Series | None = None, hydro: ~pandas.core.series.Series | None = None, other_renewables: ~pandas.core.series.Series | None = None, import_cap: ~pandas.core.series.Series | None = None, import_price: ~pandas.core.series.Series | None = None, export_cap: ~pandas.core.series.Series | None = None, export_price: ~pandas.core.series.Series | None = None, phi_fix_t: ~pandas.core.series.Series | None = None, phi_var_t: ~pandas.core.series.Series | None = None, month_of_hour: ~pandas.core.series.Series | None = None, scenario_id: int = 1, year: int = 2030, formulation_map: dict[str, str] = <factory>, cem_data: dict | None = None)[source]

Fixed-capacity designed system loaded from SDOM output snapshots.

Parameters:

  • storage_caps (dict) – Mapping {tech: {"Cap_Pch", "Cap_Pdis", "Cap_E", "eta_ch", "eta_dis", "soc_min_frac", "vom", "fom", "cost_ratio"}} for each storage technology with non-zero capacity. Capacities are in MW / MWh. fom is the fixed-O&M rate (USD/kW-yr) and cost_ratio is the share of that rate billed against the charge-side power capacity (the remainder is billed against the discharge-side power capacity), matching the CEM accounting in sdom.models.formulations_storage.storage_fixed_om_cost_expr_rule().

  • thermal_caps (dict) – Mapping {tech: {"capacity_MW", "heat_rate", "fuel_cost", "vom", "var_cost", "fom"}} for each thermal technology with non-zero capacity. var_cost = heat_rate * fuel_cost + vom. fom is the fixed-O&M rate (USD/kW-yr) aggregated across plants.

  • solar_caps (dict) – Mapping {plant_id: capacity_MW} for selected solar plants.

  • wind_caps (dict) – Mapping {plant_id: capacity_MW} for selected wind plants.

  • solar_fom (dict) – Mapping {plant_id: fom_USD_per_kW_yr} for the selected solar / wind plants. Values come from the FOM_M column of the CapSolar_*.csv / CapWind_*.csv previous-stage inputs. Carried for auditing and downstream reporting; the baseline-dispatch objective itself sources FOM from the CEM block expressions so the two paths cannot diverge.

  • wind_fom (dict) – Mapping {plant_id: fom_USD_per_kW_yr} for the selected solar / wind plants. Values come from the FOM_M column of the CapSolar_*.csv / CapWind_*.csv previous-stage inputs. Carried for auditing and downstream reporting; the baseline-dispatch objective itself sources FOM from the CEM block expressions so the two paths cannot diverge.

  • load (pandas.Series) – Hourly time-series (length 8760) indexed by hour-of-year (1..8760).

  • nuclear (pandas.Series) – Hourly time-series (length 8760) indexed by hour-of-year (1..8760).

  • hydro (pandas.Series) – Hourly time-series (length 8760) indexed by hour-of-year (1..8760).

  • other_renewables (pandas.Series) – Hourly time-series (length 8760) indexed by hour-of-year (1..8760).

  • cf_solar (pandas.DataFrame) – Hourly capacity factors with columns indexed by plant id.

  • cf_wind (pandas.DataFrame) – Hourly capacity factors with columns indexed by plant id.

  • import_cap (pandas.Series) – Hourly grid-exchange capacity and price series.

  • import_price (pandas.Series) – Hourly grid-exchange capacity and price series.

  • export_cap (pandas.Series) – Hourly grid-exchange capacity and price series.

  • export_price (pandas.Series) – Hourly grid-exchange capacity and price series.

  • phi_fix_t (pandas.Series) – Hourly fixed and variable demand-charge tariffs (USD/MW or USD/MWh).

  • phi_var_t (pandas.Series) – Hourly fixed and variable demand-charge tariffs (USD/MW or USD/MWh).

  • month_of_hour (pandas.Series) – Mapping from hour-of-year (1..8760) to calendar month (1..12) used to bill demand charges per month.

  • scenario_id (int) – Scenario / Run id resolved from the snapshot CSVs.

  • year (int) – Calendar year of the snapshot.

  • formulation_map (dict) – Mapping {component: formulation_name} resolved from defaults plus user-provided overrides.

  • cem_data (dict, optional) – CEM-shaped data dict (as returned by sdom.io_manager.load_data()) used by the baseline dispatch builder to reuse the planning-model formulations in sdom.models. None when the previous-stage inputs were not reloaded for that purpose.

storage_caps: dict[str, dict[str, float]]
thermal_caps: dict[str, dict[str, float]]
solar_caps: dict[str, float]
wind_caps: dict[str, float]
solar_fom: dict[str, float]
wind_fom: dict[str, float]
load: Series | None = None
cf_solar: DataFrame | None = None
cf_wind: DataFrame | None = None
nuclear: Series | None = None
hydro: Series | None = None
other_renewables: Series | None = None
import_cap: Series | None = None
import_price: Series | None = None
export_cap: Series | None = None
export_price: Series | None = None
phi_fix_t: Series | None = None
phi_var_t: Series | None = None
month_of_hour: Series | None = None
scenario_id: int = 1
year: int = 2030
formulation_map: dict[str, str]
cem_data: dict | None = None
__init__(storage_caps: dict[str, dict[str, float]] = <factory>, thermal_caps: dict[str, dict[str, float]] = <factory>, solar_caps: dict[str, float] = <factory>, wind_caps: dict[str, float] = <factory>, solar_fom: dict[str, float] = <factory>, wind_fom: dict[str, float] = <factory>, load: ~pandas.core.series.Series | None = None, cf_solar: ~pandas.core.frame.DataFrame | None = None, cf_wind: ~pandas.core.frame.DataFrame | None = None, nuclear: ~pandas.core.series.Series | None = None, hydro: ~pandas.core.series.Series | None = None, other_renewables: ~pandas.core.series.Series | None = None, import_cap: ~pandas.core.series.Series | None = None, import_price: ~pandas.core.series.Series | None = None, export_cap: ~pandas.core.series.Series | None = None, export_price: ~pandas.core.series.Series | None = None, phi_fix_t: ~pandas.core.series.Series | None = None, phi_var_t: ~pandas.core.series.Series | None = None, month_of_hour: ~pandas.core.series.Series | None = None, scenario_id: int = 1, year: int = 2030, formulation_map: dict[str, str] = <factory>, cem_data: dict | None = None) → None
class sdom.resiliency.BaselineState(soc_trajectory: ~pandas.core.frame.DataFrame | None = None, solver_status: str | None = None, objective_value: float | None = None, metadata: dict[str, ~typing.Any] = <factory>)[source]

Placeholder container for baseline-dispatch outputs (Phase 2).

Parameters:

  • soc_trajectory (pandas.DataFrame, optional) – Hourly state-of-charge per storage technology (hour x tech).

  • solver_status (str, optional) – Solver termination status from the baseline run.

  • objective_value (float, optional) – Baseline objective value (USD).

  • metadata (dict, optional) – Free-form solver / run metadata.

soc_trajectory: DataFrame | None = None
solver_status: str | None = None
objective_value: float | None = None
metadata: dict[str, Any]
__init__(soc_trajectory: ~pandas.core.frame.DataFrame | None = None, solver_status: str | None = None, objective_value: float | None = None, metadata: dict[str, ~typing.Any] = <factory>) → None

Outage Specification

class sdom.resiliency.OutageSpec(duration_hours: int, recovery_hours: int | dict[str, int], outaged_assets: dict[str, str | ~typing.Iterable], derating_factors: dict[tuple[str, str], float] = <factory>, min_soc_recovery: dict[str, float] | None = None, per_asset_durations: dict[tuple[str, str], int] = <factory>)[source]

Specification of an outage / de-rating scenario.

Parameters:

  • duration_hours (int) – Outage duration applied to all listed assets unless overridden by per_asset_durations.

  • recovery_hours (int or dict) – Hours allowed after the outage for storage devices to recover the SOC target. A single int broadcasts to every storage technology; a dict {tech: hours} gives per-tech values.

  • outaged_assets (dict) – Mapping {component: asset_selector} where component is one of VALID_COMPONENTS and asset_selector is either the string "all" or an iterable of asset/plant IDs. The must-run components in MUST_RUN_COMPONENTS only accept "all" in this iteration; iterables raise NotImplementedError from validate().

  • derating_factors (dict, optional) – Mapping {(component, asset_id): factor} with each factor in [0, 1]. Missing assets default to 0 (full outage) when listed in outaged_assets, or 1 (no outage) otherwise.

  • min_soc_recovery (dict, optional) – Per-tech target SOC fraction at the end of the recovery window. None (default) -> baseline SOC at the end of the recovery window divided by Cap_E per tech.

  • per_asset_durations (dict, optional) – Optional per-asset overrides {(component, asset_id): hours} for duration_hours.

Raises:

ValueError – If a derating factor is outside [0, 1], a component name is not in VALID_COMPONENTS, or a string asset selector is not "all".

duration_hours: int
recovery_hours: int | dict[str, int]
outaged_assets: dict[str, str | Iterable]
derating_factors: dict[tuple[str, str], float]
min_soc_recovery: dict[str, float] | None = None
per_asset_durations: dict[tuple[str, str], int]
__post_init__() → None[source]
validate(designed_system) → None[source]

Validate the specification against designed_system.

Parameters:

designed_system (DesignedSystem) – Source of truth for asset universes (storage techs, thermal plants, VRE plants).

Raises:

  • ValueError – If an asset id in outaged_assets does not match the designed system, or if recovery_hours is a dict missing techs present in the designed system.

  • NotImplementedError – If a must-run component is given a per-asset iterable.

resolve_recovery_hours(designed_system) → dict[str, int][source]

Return {tech: recovery_hours} for every storage tech.

Raises:

ValueError – If recovery_hours is a dict that does not cover every storage tech in designed_system.

resolve_duration(component: str, asset_id: str) → int[source]

Return outage duration (hours) for a given (component, asset_id).

resolve_derating(component: str, asset_id: str) → float[source]

Return the derating multiplier rho for (component, asset_id).

Returns 1.0 when the asset is not selected for outage. Returns the user-provided factor (or the default 0.0) when it is.

resolve_outaged_asset_ids(component: str, designed_system) → list[str][source]

Return the concrete list of outaged asset ids for component.

For must-run components and imports the returned list is the canonical universe (["all"] and ["grid"] respectively).

resolve_min_soc_recovery(baseline_results, designed_system, *, recovery_end_hour: dict[str, int]) → dict[str, float][source]

Return {tech: SOC_target_fraction} at the end of recovery.

If min_soc_recovery is None, fractions default to SOC_baseline[tech, recovery_end_hour[tech]] / Cap_E[tech]. Per-tech entries in min_soc_recovery override the baseline default.

Parameters:

  • baseline_results (BaselineDispatchResults) – Provides soc_trajectory (DataFrame indexed by hour with one column per tech, MWh).

  • designed_system (DesignedSystem) – Provides storage_caps[tech]["Cap_E"].

  • recovery_end_hour (dict) – Mapping {tech: hour_of_year} of the recovery end-hour (already clipped to the baseline horizon).

sdom.resiliency.VALID_COMPONENTS

Built-in immutable sequence.

If no argument is given, the constructor returns an empty tuple. If iterable is specified the tuple is initialized from iterable’s items.

If the argument is a tuple, the return value is the same object.

sdom.resiliency.MUST_RUN_COMPONENTS

Built-in immutable sequence.

If no argument is given, the constructor returns an empty tuple. If iterable is specified the tuple is initialized from iterable’s items.

If the argument is a tuple, the return value is the same object.

Baseline Dispatch (Problem B)

sdom.resiliency.build_baseline_dispatch(designed_system, *, n_hours=8760, min_soc_per_tech=None, curtailment_penalty=0.0, add_demand_charges=True, model_name='SDOM_BaselineDispatch', profile=False)[source]

Build the fixed-capacity annual baseline dispatch model.

The model is built by calling the planning-model constructor sdom.optimization_main._initialize_model_copperplate() on designed_system.cem_data and then:

  1. Fixing every capacity decision variable (model.thermal.plant_installed_capacity, model.storage.{Pcha,Pdis,Ecap}, model.{pv,wind}.capacity_fraction) to the values stored in designed_system.

  2. Applying an optional operational SOC floor.

  3. Layering monthly fixed and variable demand charges on top of the existing model.imports.variable via add_demand_charges_to_existing_imports().

  4. Replacing the CEM planning objective with an operational-only objective (variable thermal + variable storage + imports/exports + demand charges + optional curtailment penalty).

Parameters:

  • designed_system (DesignedSystem) – Output of sdom.resiliency.load_designed_system(). Must carry a populated cem_data attribute.

  • n_hours (int, optional) – Number of hours to simulate. Default 8760.

  • min_soc_per_tech (dict, optional) – {tech: fraction in [0, 1]} enforcing a SOC floor as a fraction of energy capacity. Missing techs default to whatever is stored on designed_system.storage_caps[s]['soc_min_frac'] (or 0).

  • curtailment_penalty (float, optional) – Per-MWh penalty applied to total VRE curtailment. Default 0.0.

  • add_demand_charges (bool, optional) – When True (default) and designed_system carries phi_fix_t, phi_var_t and month_of_hour, monthly demand charges are added to the objective.

  • model_name (str, optional) – Pyomo model name. Default "SDOM_BaselineDispatch".

  • profile (bool, optional) – Instrument every build step with a ModelInitProfiler and print a summary. Default False.

Return type:

pyomo.environ.ConcreteModel

Raises:
sdom.resiliency.run_baseline_dispatch(model, *, solver='highs', solver_options=None, tee=False, profile=False)[source]

Solve the baseline dispatch and collect per-hour trajectories.

Parameters:

  • model (pyomo.environ.ConcreteModel) – Model returned by build_baseline_dispatch().

  • solver (str, optional) – Pyomo solver name. "highs" first tries appsi_highs then falls back to highs. Default "highs".

  • solver_options (dict, optional) – Extra options passed to solver.solve(..., options=...).

  • tee (bool, optional) – Stream solver output to the console. Default False.

  • profile (bool, optional) – Instrument solve + extraction with a profiler. Default False.

Returns:

Hourly dispatch trajectories plus solver metadata.

Return type:

BaselineDispatchResults

Raises:

AttributeError – If model was not produced by build_baseline_dispatch().

class sdom.resiliency.BaselineDispatchResults(soc_trajectory: ~pandas.core.frame.DataFrame | None = None, pcha_trajectory: ~pandas.core.frame.DataFrame | None = None, pdis_trajectory: ~pandas.core.frame.DataFrame | None = None, pthermal_trajectory: ~pandas.core.frame.DataFrame | None = None, psolar_trajectory: ~pandas.core.frame.DataFrame | None = None, pwind_trajectory: ~pandas.core.frame.DataFrame | None = None, pimp: ~pandas.core.series.Series | None = None, pexp: ~pandas.core.series.Series | None = None, nuclear: ~pandas.core.series.Series | None = None, hydro: ~pandas.core.series.Series | None = None, other_renewables: ~pandas.core.series.Series | None = None, load: ~pandas.core.series.Series | None = None, month_of_hour: ~pandas.core.series.Series | None = None, objective_value: float | None = None, solver_status: str | None = None, metadata: dict[str, ~typing.Any] = <factory>, cost_breakdown: dict[str, float] = <factory>)[source]

Trajectories and metadata produced by run_baseline_dispatch().

Parameters:

  • soc_trajectory (pandas.DataFrame) – Hourly state-of-charge per storage technology, indexed by hour and with one column per tech (MWh).

  • pcha_trajectory (pandas.DataFrame) – Hourly charge / discharge per storage tech (MW).

  • pdis_trajectory (pandas.DataFrame) – Hourly charge / discharge per storage tech (MW).

  • pthermal_trajectory (pandas.DataFrame) – Hourly thermal dispatch per balancing-unit Plant_id (MW). Empty DataFrame when no thermal units survive the snapshot filter.

  • psolar_trajectory (pandas.DataFrame) – Hourly dispatched solar / wind power per plant id (MW).

  • pwind_trajectory (pandas.DataFrame) – Hourly dispatched solar / wind power per plant id (MW).

  • pimp (pandas.Series) – Hourly imports / exports (MW).

  • pexp (pandas.Series) – Hourly imports / exports (MW).

  • nuclear (pandas.Series) – Hourly time-series parameters echoed from the input system (MW).

  • hydro (pandas.Series) – Hourly time-series parameters echoed from the input system (MW).

  • other_renewables (pandas.Series) – Hourly time-series parameters echoed from the input system (MW).

  • load (pandas.Series) – Hourly time-series parameters echoed from the input system (MW).

  • month_of_hour (pandas.Series) – Hour -> month mapping used by the demand-charge billing.

  • objective_value (float) – Operational objective value (USD).

  • solver_status (str) – Solver termination condition (e.g. "optimal").

  • metadata (dict, optional) – Free-form solver / run metadata.

  • cost_breakdown (dict, optional) – Per-component USD totals reconciling to objective_value. Keys: thermal_var_USD, storage_var_USD, imports_USD, exports_USD (positive; objective contribution is -exports_USD), demand_charges_USD, curtailment_USD, fom_USD, total_USD. Empty dict when the model carries no component metadata.

soc_trajectory: DataFrame | None = None
pcha_trajectory: DataFrame | None = None
pdis_trajectory: DataFrame | None = None
pthermal_trajectory: DataFrame | None = None
psolar_trajectory: DataFrame | None = None
pwind_trajectory: DataFrame | None = None
pimp: Series | None = None
pexp: Series | None = None
nuclear: Series | None = None
hydro: Series | None = None
other_renewables: Series | None = None
load: Series | None = None
month_of_hour: Series | None = None
objective_value: float | None = None
solver_status: str | None = None
metadata: dict[str, Any]
cost_breakdown: dict[str, float]
__init__(soc_trajectory: ~pandas.core.frame.DataFrame | None = None, pcha_trajectory: ~pandas.core.frame.DataFrame | None = None, pdis_trajectory: ~pandas.core.frame.DataFrame | None = None, pthermal_trajectory: ~pandas.core.frame.DataFrame | None = None, psolar_trajectory: ~pandas.core.frame.DataFrame | None = None, pwind_trajectory: ~pandas.core.frame.DataFrame | None = None, pimp: ~pandas.core.series.Series | None = None, pexp: ~pandas.core.series.Series | None = None, nuclear: ~pandas.core.series.Series | None = None, hydro: ~pandas.core.series.Series | None = None, other_renewables: ~pandas.core.series.Series | None = None, load: ~pandas.core.series.Series | None = None, month_of_hour: ~pandas.core.series.Series | None = None, objective_value: float | None = None, solver_status: str | None = None, metadata: dict[str, ~typing.Any] = <factory>, cost_breakdown: dict[str, float] = <factory>) → None

Outage Dispatch (Problem O)

sdom.resiliency.build_outage_dispatch(baseline_results, *, start_hour, outage_spec, designed_system=None, slack_penalty=10000.0, curtailment_penalty=0.0, soc_slack_penalty=1000.0, min_soc_per_tech=None, n_hours=8760, critical_load_MW=None, model_name='SDOM_OutageDispatch', profile=False)[source]

Build the per-hour outage economic-dispatch Pyomo LP anchored at start_hour.

Parameters:

  • baseline_results (BaselineDispatchResults) – Output of sdom.resiliency.run_baseline_dispatch(). Used to seed the initial SOC and (optionally) the recovery SOC target. May carry the originating DesignedSystem under baseline_results.metadata["designed_system"].

  • start_hour (int) – Anchor hour h. The outage horizon is [h, h + duration + max_recovery - 1], clipped to [1, n_hours].

  • outage_spec (OutageSpec) – Outage / de-rating specification.

  • designed_system (DesignedSystem, optional) – Source of truth for capacities and time series. If None, taken from baseline_results.metadata["designed_system"].

  • slack_penalty (float, optional) – Penalty \(\pi^{slack}\) (USD/MWh) on unserved-energy slack. Default 10_000.0.

  • curtailment_penalty (float, optional) – Penalty on curtailed VRE energy (USD/MWh). Default 0.0.

  • soc_slack_penalty (float, optional) – Penalty \(\pi^{soc}\) (USD/MWh) on the per-storage-tech slack variable that relaxes the SOC recovery-target constraint (see notes). Default 1_000.0. The operational SOC floor remains a hard bound; only the end-of-recovery target is relaxed.

  • min_soc_per_tech (dict, optional) – Operational SOC floor per storage tech (fraction of Cap_E). Same semantics as build_baseline_dispatch().

  • n_hours (int, optional) – Length of the baseline horizon used for end-of-year clipping. Default 8760.

  • critical_load_MW (float, optional) – Constant critical load (MW) used in place of designed_system.load[t] for every hour t in the outage sub-horizon [start_hour, start_hour + duration_hours - 1] (clipped to the LP end hour). Recovery-window hours continue to use the original D_t. None (default) preserves the original behaviour: the hourly load series is used everywhere. Must be non-negative.

  • model_name (str, optional) – Pyomo model name. Default "SDOM_OutageDispatch".

  • profile (bool, optional) – When True, instrument the build stages with a ModelInitProfiler, attach it as model.profiler and print a summary table. Default False. Note: enabling this from inside a parallel ProcessPoolExecutor will produce one summary per worker on each spawned process, which is rarely useful. Prefer profile only in serial runs.

Returns:

A Pyomo LP exposing model.h (hour set), model.u (slack, NonNegativeReals), the standard dispatch sub-blocks (storage, thermal, solar, wind, imports, exports) without demand-charge variables, model.fom_cost_expr (constant fixed-O&M cost prorated to the outage horizon), and an objective that minimises operational cost plus slack and curtailment penalties plus the prorated FOM constant.

Return type:

pyomo.environ.ConcreteModel

Raises:

  • ValueError – If designed_system is not provided and cannot be recovered from baseline_results.metadata, or if validation of outage_spec against designed_system fails.

  • TypeError – If baseline_results or designed_system is the wrong type.

Notes

Initial SOC is seeded by setting the mutable parameter model.storage.SOC_init[s] from baseline_results.soc_trajectory at start_hour. SOC_init represents the SOC at the start of the outage horizon (i.e., the boundary value \(SOC_{s,h-1}\) that the dynamics equation at \(t = h\) reads as its prior state). The SOC dynamics constraint therefore covers every hour in \(\mathcal{T}^{out}_h\), including the anchor hour start_hour, so that the charge and discharge variables at the anchor hour appear in a balance equation. Earlier versions fixed SOC[s, start_hour] via pyomo.environ.Var.fix() and skipped the dynamics equation at the anchor; under that formulation Pcha[s, start_hour] and Pdis[s, start_hour] for surviving (non-outaged) storage techs were unconstrained by any SOC balance.

When critical_load_MW is provided, the load parameter is overridden only over the outage sub-horizon \(\mathcal{T}^{out}_h = \{h, \ldots, h + \Delta^{out} - 1\}\); hours in the recovery sub-horizon retain the original D_t so that storage replenishment toward the end-of-recovery target reflects realistic post-outage operations.

Imports Formulation with Demand Charges

The opt-in ImportsWithDemandChargesFormulation block builder used by the baseline LP. Pure linear program; not registered in io_manager.py.

sdom.resiliency.add_imports_with_demand_charges(model, *, import_cap: Series, import_price: Series, phi_fix_t: Series, phi_var_t: Series, month_of_hour: Series, block_name: str = 'imports')[source]

Attach the ImportsWithDemandChargesFormulation block to model.

The block adds hourly import variables Pimp[t], monthly fixed and variable demand-charge variables D_fix[m] / D_var[m], the capacity bound, the demand-charge linking inequalities, and an additive cost expression total_cost_expr suitable for inclusion in any objective.

Parameters:

  • model (pyomo.environ.ConcreteModel) – Host model. Must already define model.h (a Pyomo Set of hour indices).

  • import_cap (pandas.Series) – Hourly import capacity \(\overline{P}^{imp}_{t}\) (MW), indexed by hour matching model.h.

  • import_price (pandas.Series) – Hourly import energy price \(c^{imp}_{t}\) (USD/MWh).

  • phi_fix_t (pandas.Series) – Hourly fixed demand-charge tariff \(\phi^{fix}_{t}\) (USD/MW). Must be constant within each calendar month.

  • phi_var_t (pandas.Series) – Hourly variable (time-of-use) demand-charge tariff \(\phi^{var}_{t}\) (USD/MW).

  • month_of_hour (pandas.Series) – Mapping hour -> month integer.

  • block_name (str, optional) – Name of the sub-block attached to model. Default "imports".

Returns:

The block that was attached to model.

Return type:

pyomo.environ.Block

Raises:

AttributeError – If model.h is not present.

Notes

The contributed cost is

\[Z_{imp,dc} = \sum_t c^{imp}_t\, p^{imp}_t + \sum_m \left( D^{fix}_m + D^{var}_m \right),\]

with linking constraints \(D^{k}_{m} \ge \phi^{k}_{t}\, p^{imp}_{t}\) for all \(t \in \mathcal{T}_m\) and \(k \in \{fix, var\}\).

Examples

>>> import pyomo.environ as pyo
>>> import pandas as pd
>>> m = pyo.ConcreteModel()
>>> m.h = pyo.RangeSet(1, 24)
>>> idx = range(1, 25)
>>> add_imports_with_demand_charges(  
...     m,
...     import_cap=pd.Series(100.0, index=idx),
...     import_price=pd.Series(1.0, index=idx),
...     phi_fix_t=pd.Series(50.0, index=idx),
...     phi_var_t=pd.Series(2.0, index=idx),
...     month_of_hour=pd.Series(1, index=idx),
... )

Parallel Runner

sdom.resiliency.run_resiliency_evaluation(baseline_results, *, outage_spec, designed_system=None, hours=None, slack_penalty=10000.0, curtailment_penalty=0.0, soc_slack_penalty=1000.0, min_soc_per_tech=None, n_hours=8760, n_workers=None, solver='highs', solver_options=None, critical_load_MW=None, profile_outages=False)[source]

Run the per-hour outage evaluation in parallel and aggregate metrics.

For every anchor hour h in hours the runner builds the short-horizon outage LP via build_outage_dispatch(), solves it with HiGHS, and records EUE, USE_hours, max_unserved_MW and the objective value. Per-hour problems are independent and are solved in parallel via concurrent.futures.ProcessPoolExecutor when n_workers > 1; n_workers == 1 (or a single-hour run) falls back to an in-process serial loop.

Parameters:

  • baseline_results (BaselineDispatchResults) – Output of run_baseline_dispatch(). Must carry an SOC trajectory; if its metadata does not contain a "designed_system" entry, designed_system must be supplied.

  • outage_spec (OutageSpec) – Outage / de-rating specification, broadcast to every anchor hour.

  • designed_system (DesignedSystem, optional) – Source of truth for capacities and time series. Takes precedence over baseline_results.metadata['designed_system'] when both are provided.

  • hours (iterable of int, optional) – Anchor hours to evaluate. None (default) evaluates every hour 1..n_hours. Duplicates are removed; the result is sorted.

  • slack_penalty (float, optional) – Penalty (USD/MWh) applied to slack u[t]. Default 10_000.

  • curtailment_penalty (float, optional) – Penalty applied to curtailed VRE energy (USD/MWh). Default 0.

  • soc_slack_penalty (float, optional) – Penalty (USD/MWh) applied to the per-tech SOC recovery-target slack variable in the outage LP. Forwarded to build_outage_dispatch(). Default 1_000.

  • min_soc_per_tech (dict, optional) – Operational SOC floor per storage tech (fraction of Cap_E); forwarded to the outage builder.

  • n_hours (int, optional) – Length of the baseline horizon used for end-of-year clipping. Default 8760.

  • n_workers (int, optional) – Worker pool size. None -> max(1, os.cpu_count() - 1); 1 forces serial mode. Always clamped to len(hours).

  • solver (str, optional) – Pyomo solver name. "highs" first tries appsi_highs. Default "highs".

  • solver_options (dict, optional) – Solver options forwarded to solver.solve(..., options=...).

  • critical_load_MW (float, optional) – Constant critical load (MW) used in place of the hourly load series during the outage sub-horizon of each per-hour LP. Forwarded to build_outage_dispatch(). None (default) preserves the original behaviour. Must be non-negative.

  • profile_outages (bool, optional) – When True and n_workers == 1, profile every per-hour outage build via ModelInitProfiler. Ignored (with a warning) when n_workers > 1 because each worker would emit its own summary on a separate process. Default False.

Returns:

Per-hour records (sorted by start_hour) plus run metadata including n_workers_used and a reference to outage_spec.

Return type:

ResiliencyResults

Raises:

Notes

Failure isolation: if any worker raises an exception, that hour’s record is marked solver_status="error" and its error_message field carries the formatted traceback. Other hours continue normally.

Results Container

class sdom.resiliency.ResiliencyResults(per_hour: ~pandas.core.frame.DataFrame, metadata: dict[str, ~typing.Any] = <factory>)[source]

Per-hour outage outcomes (lightweight, Phase 5).

Aggregate metrics (LOLP, LOLE, percentiles) and plotting are added in Phase 6.

Parameters:

  • per_hour (pandas.DataFrame) – Indexed by hour (anchor start_hour). Columns include ["EUE", "USE_hours", "max_unserved_MW", "objective_value", "solver_status", "solve_time_s", "truncated", "error_message"].

  • metadata (dict) – Free-form run metadata. Conventionally includes {"n_workers_used", "outage_spec", "n_hours", "solver"}.

per_hour: DataFrame
metadata: dict[str, Any]
to_dataframe() → DataFrame[source]

Return per_hour with the index promoted to a hour column.

Returns:

A copy of per_hour with hour as a regular column, sorted by hour.

Return type:

pandas.DataFrame

eue_total() → float[source]

Return the sum of per-hour expected unserved energy (MWh).

Return type:

float

metrics(*, level: str = 'aggregate')[source]

Aggregate or per-hour resiliency metrics.

Parameters:

level ({"aggregate", "per_hour"}, optional) – "aggregate" (default) returns a dict of scalar metrics computed over the evaluated hours (errored hours excluded). "per_hour" returns a copy of per_hour with hour promoted to a column.

Return type:

dict or pandas.DataFrame

Raises:

ValueError – If level is not one of the supported values.

Notes

Aggregate metrics exclude rows with solver_status == "error"; the count of excluded rows is reported as n_errors.

lolp() → float[source]

Return the loss-of-load probability across evaluated hours.

Return type:

float

lole() → float[source]

Return the loss-of-load expectation (mean USE hours per scenario).

Return type:

float

eue(*, p: float | None = None) → float[source]

Return the mean EUE or an empirical percentile of EUE.

Parameters:

p (float, optional) – Quantile in (0, 1). Default None returns the mean EUE.

Return type:

float

Raises:

ValueError – If p is provided and not in (0, 1).

save(path: str | Path | None = None) → Path[source]

Persist per-hour records and aggregate metrics to disk.

Parameters:

path (str or pathlib.Path, optional) – Output directory. Default: ./results_resiliency/ relative to the current working directory. The directory is created if it does not exist.

Returns:

The directory the artifacts were written to.

Return type:

pathlib.Path

Raises:

ImportError – If no Parquet engine (pyarrow or fastparquet) is available.

Notes

Writes two files to path:

  • per_hour.parquet - the per-hour DataFrame.

  • summary.json - aggregate metrics + JSON-safe metadata.

classmethod load(path: str | Path) → ResiliencyResults[source]

Load a previously-saved ResiliencyResults from path.

Parameters:

path (str or pathlib.Path) – Directory that previously received save().

Return type:

ResiliencyResults

Raises:

FileNotFoundError – If per_hour.parquet or summary.json is missing.

__init__(per_hour: ~pandas.core.frame.DataFrame, metadata: dict[str, ~typing.Any] = <factory>) → None

Plotting

sdom.resiliency.plot_metric_distribution(results: ResiliencyResults, *, metric: str = 'EUE', kind: str = 'hist', ax: Axes | None = None, **plot_kwargs: Any) → Axes[source]

Plot the empirical distribution of a per-hour metric.

Parameters:

  • results (ResiliencyResults) – Container produced by sdom.resiliency.run_resiliency_evaluation().

  • metric (str, optional) – Numeric column of results.per_hour to plot. Default "EUE".

  • kind ({"hist", "ecdf", "exceedance"}, optional) –

    Plot style. Default "hist".

    • "hist" - histogram of metric values.

    • "ecdf" - empirical CDF, monotonically non-decreasing in [0, 1].

    • "exceedance" - exceedance curve 1 - ECDF, monotonically non-increasing.

  • ax (matplotlib.axes.Axes, optional) – Existing axes to draw on. A new figure/axes is created when None.

  • **plot_kwargs – Forwarded to the underlying matplotlib call (ax.hist for kind="hist"; ax.plot otherwise).

Return type:

matplotlib.axes.Axes

Raises:

  • ImportError – If matplotlib is not importable.

  • ValueError – If kind is not a supported value or metric is not a numeric column of results.per_hour.