How to Define and Apply Contingencies

To follow along, you can download this tutorial as a Julia script (.jl) or Jupyter notebook (.ipynb).

This guide shows how to compute post-contingency PTDF rows with VirtualMODF — the lazy Multiple Outage Distribution Factor matrix. You will attach outages to a system, let them auto-register, query monitored arcs under a contingency, and build manual modifications when you need full control.

Prerequisites

• PowerNetworkMatrices.jl and PowerSystems.jl installed
• A power system model

using PowerNetworkMatrices
import PowerSystems
import PowerSystemCaseBuilder

sys = PowerSystemCaseBuilder.build_system(
    PowerSystemCaseBuilder.PSITestSystems,
    "c_sys5",
);

[ Info: Loaded time series from storage file existing=/home/runner/.julia/packages/PowerSystemCaseBuilder/zuCHR/data/serialized_system/8584b9e729c8aa68ee5405660c6258cde1f67ed3b68f114823707912e9a0d16c/c_sys5_time_series_storage.h5 new=/tmp/jl_ELSW9f compression=InfrastructureSystems.CompressionSettings(false, InfrastructureSystems.CompressionTypesModule.CompressionTypes.DEFLATE = 1, 3, true)

Attach Outages to the System

Contingencies are defined as PowerSystems.Outage supplemental attributes on the components they trip. When a contingency only needs to exist, use a PowerSystems.FixedForcedOutage with outage_status = 1.0 (outaged) and attach it to each branch:

for branch in PowerSystems.get_components(PowerSystems.ACTransmission, sys)
    outage = PowerSystems.FixedForcedOutage(; outage_status = 1.0)
    PowerSystems.add_supplemental_attribute!(sys, branch, outage)
end

Build the VirtualMODF

Registration is automatic: the VirtualMODF constructor scans the system for PowerSystems.Outage attributes and resolves each into a ContingencySpec.

vmodf = VirtualMODF(sys)

VirtualMODF{Tuple{Vector{Tuple{Int64, Int64}}, Vector{Int64}}, Tuple{Dict{Tuple{Int64, Int64}, Int64}, Dict{Int64, Int64}}, PowerNetworkMatrices.KLUWrapper.KLULinSolveCache{Float64, Int64}}
    Dimension 1, [(1, 2), (1, 5), (3, 4), (2, 3), (4, 5), (1, 4)]
    Dimension 2, [1, 2, 3, 4, 5]
And data with size (6, 5):
VirtualMODF with 6 registered contingencies

Inspect what was registered with get_registered_contingencies. It returns a Dict{UUID, ContingencySpec} keyed by the source outage's UUID:

contingencies = get_registered_contingencies(vmodf)

Dict{Base.UUID, ContingencySpec} with 6 entries:
  UUID("2429f852-5b99-449a… => ContingencySpec(UUID("2429f852-5b99-449a-a7a0-ba…
  UUID("85c374dc-f3d1-4b62… => ContingencySpec(UUID("85c374dc-f3d1-4b62-a049-59…
  UUID("de8f5d2a-5f03-4fdd… => ContingencySpec(UUID("de8f5d2a-5f03-4fdd-ad86-91…
  UUID("fee1ff1e-c004-40ea… => ContingencySpec(UUID("fee1ff1e-c004-40ea-9a12-d9…
  UUID("ecc45dfa-8e93-4b27… => ContingencySpec(UUID("ecc45dfa-8e93-4b27-a5ce-9a…
  UUID("c24ce884-c7f1-4a84… => ContingencySpec(UUID("c24ce884-c7f1-4a84-b326-51…

Query a Monitored Arc Under a Contingency

Index the matrix as vmodf[monitored_arc, spec]. The monitored arc is an arc tuple (from, to) (or its integer index); the returned value is the full post-contingency PTDF row for that arc — one sensitivity per bus.

Pick a monitored arc and outage a different arc — monitoring an element that the contingency itself outages is undefined and raises. get_arc_axis lists the arcs available to pick from, and the convenience NetworkModification constructor builds the spec straight from one of them:

arcs = get_arc_axis(vmodf)
monitored_arc = arcs[1]
ctg = NetworkModification(vmodf, arcs[2])

NetworkModification("outage_1_5", (ArcModification(2, -156.25, -15.470297f0 + 154.68735f0im, 15.470297f0 - 154.70297f0im, 15.470297f0 - 154.70297f0im, -15.470297f0 + 154.68735f0im),), (), false)

The returned row carries one post-contingency sensitivity per bus:

row = vmodf[monitored_arc, ctg]

5-element Vector{Float64}:
  0.3070707007235925
 -0.40909091870631104
 -0.30000000787526937
  0.0
  1.5433482929455187e-16

The second index accepts three equivalent forms — the NetworkModification used above, a ContingencySpec from the registered set, or the original PowerSystems.Outage (by its registered UUID). All resolve to the same NetworkModification and share the cached Woodbury factors, so repeated queries for one contingency across different monitored arcs reuse work.

A ContingencySpec pairs a NetworkModification with the UUID of the outage it came from, so the registered spec matching ctg is the one whose modification field compares equal to it:

spec = first(s for s in values(contingencies) if s.modification == ctg)
vmodf[monitored_arc, spec] == vmodf[monitored_arc, spec.modification] == row

true

Indexing by the PowerSystems.Outage attribute itself gives the same row again. Fetch the outage from the branch that ctg trips:

outaged_branch = first(
    branch for branch in PowerSystems.get_components(PowerSystems.ACTransmission, sys)
    if PowerSystems.get_number(PowerSystems.get_from(PowerSystems.get_arc(branch))) ==
    arcs[2][1] &&
    PowerSystems.get_number(PowerSystems.get_to(PowerSystems.get_arc(branch))) ==
    arcs[2][2]
)
outage =
    first(PowerSystems.get_supplemental_attributes(PowerSystems.Outage, outaged_branch))
vmodf[monitored_arc, outage] == row

true

The modification type model

Under the convenience constructor sit a few value types, layered from the smallest unit up to the solver-ready form. It helps to know them before dropping to the manual path:

NetworkModification is the canonical representation: once built it holds no reference to the source System and serves as the cache key inside VirtualMODF (its label is excluded from equality, so two physically identical modifications compare equal regardless of name).

Build a Modification Manually

The convenience constructor used above (NetworkModification(matrix, arc), or NetworkModification(matrix, branch)) is the simplest path — it looks up the arc's susceptance and populates the deltas for you. When you want full control, assemble the low-level building blocks instead. An ArcModification is a susceptance change on one arc (delta_b negative for an outage); a ShuntModification is an admittance change on one bus. Both are indexed by their integer position in the matrix, which get_arc_lookup and get_bus_lookup provide:

arc_index = get_arc_lookup(vmodf)[(1, 4)]
arc_mod = ArcModification(arc_index, -5.0)   # a partial Δb on a direct arc

ArcModification(6, -5.0, 0.0f0 + 0.0f0im, 0.0f0 + 0.0f0im, 0.0f0 + 0.0f0im, 0.0f0 + 0.0f0im)

A ShuntModification is built the same way, from the bus position that get_bus_lookup returns:

bus_index = get_bus_lookup(vmodf)[3]
shunt_mod = ShuntModification(bus_index, ComplexF32(-0.1im))

ShuntModification(3, -0.0f0 - 0.1f0im)

Combining the two gives one NetworkModification — the constructor takes a label, the arc changes, the shunt changes, and whether the modification islands the network:

custom = NetworkModification("arc_and_shunt", [arc_mod], [shunt_mod], false)
vmodf[monitored_arc, custom]

5-element Vector{Float64}:
  0.20773367562237602
 -0.46773741395863283
 -0.3430074378450955
  0.0
  0.17090554581590875

Prefer the convenience constructors over hand-built ArcModification values: they compute physically consistent delta_b and Pi-model deltas from the network data, which is otherwise your responsibility to get right. Compare the -5.0 above with the delta_b that NetworkModification derives for a full outage of the same arc:

only(NetworkModification(vmodf, (1, 4)).arc_modifications).delta_b

-32.894737243652344

One-Shot Post-Modification Rows from a VirtualPTDF

If you already hold a VirtualPTDF and want a single post-modification row without registering contingencies, use get_post_modification_ptdf_row. It applies a NetworkModification through the same Woodbury correction, on an arc taken from get_arc_axis:

vptdf = VirtualPTDF(sys)
varcs = get_arc_axis(vptdf)
mod = NetworkModification(vptdf, varcs[2])
row_oneshot = get_post_modification_ptdf_row(vptdf, varcs[1], mod)

5-element Vector{Float64}:
  0.3070707007235925
 -0.40909091870631104
 -0.30000000787526937
  0.0
  1.5433482929455187e-16

Indexing is the equivalent form — it returns the same row:

isapprox(vptdf[varcs[1], mod], row_oneshot)

true

This function does no caching — each call recomputes. When querying many monitored arcs for the same modification, precompute once with compute_woodbury_factors and reuse via apply_woodbury_correction.

Contingencies and Network Reduction

If you build the VirtualMODF with network_reductions, any branch that a contingency outages or monitors must survive every reduction step. Outage and monitored-component buses are auto-protected from reduction. Declare monitored branches on the outage so their buses are kept:

monitored_line = PowerSystems.get_component(PowerSystems.ACTransmission, sys, "2")
PowerSystems.set_monitored_components!(outage, [monitored_line])

Set{Base.UUID} with 1 element:
  Base.UUID("46016b7c-9676-4d1e-8ebe-7b6ec640ebcb")

The reduced VirtualMODF then keeps every arc the registered contingencies need. Reductions are passed as a vector of NetworkReduction specs, and get_arc_axis reports what survived. c_sys5 has no radial buses to remove, so here RadialReduction leaves the arc axis untouched; on a network that does reduce, the protected arcs are the ones guaranteed to survive:

vmodf_reduced = VirtualMODF(sys; network_reductions = NetworkReduction[RadialReduction()])
get_arc_axis(vmodf_reduced)

6-element Vector{Tuple{Int64, Int64}}:
 (1, 2)
 (1, 5)
 (3, 4)
 (2, 3)
 (4, 5)
 (1, 4)

Querying a monitored arc that was reduced away raises a clear error rather than silently returning the base row.

See Also

• Public API Reference — full docstrings for NetworkModification, ContingencySpec, compute_woodbury_factors, and the rest
• Flowgate Methodology — the Woodbury post-contingency theory