Per-type service models and merged sparse/dense reserve containers

docs/Project.toml

@@ -8,6 +8,9 @@ InfrastructureSystems = "2cd47ed4-ca9b-11e9-27f2-ab636a7671f1"
 Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306"
 PrettyTables = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d"
 
+[sources]
+InfrastructureSystems = {url = "https://github.com/Sienna-Platform/InfrastructureSystems.jl", rev = "IS4"}
+
 [compat]
 Documenter = "^1.0"
 julia = "^1.10"

docs/make.jl

@@ -36,6 +36,11 @@ makedocs(
     pages = Any[p for p in pages],
     draft = false,
     plugins = [links],
+    # Pre-existing broken `@ref` links (also red on `main`): the explanation/reference pages
+    # reference symbols that were moved to PowerOperationsModels in the IOM/POM split (#104).
+    # Downgrade only the cross-reference failures to warnings so the site builds; every other
+    # category (missing docstrings, doctests, ...) still fails the build.
+    warnonly = [:cross_references],
 )
 
 deploydocs(

src/InfrastructureOptimizationModels.jl

@@ -297,7 +297,6 @@ export get_parameter_array
 export get_network_reduction
 export get_multiplier_array
 export get_parameter_column_refs
-export get_service_name
 export get_default_time_series_type
 export add_expression_container!
 

src/common_models/add_constraint_dual.jl

@@ -30,40 +30,32 @@ function add_constraint_dual!(
 end
 
 # Service model
+#
+# Services of the same type share merged constraint containers keyed by
+# `(constraint_type, service_type)` with empty meta, so the dual mirrors the existing
+# constraint container exactly (as the device/network paths do) rather than building a
+# per-service `[service_name]` axis. The `haskey` check keeps dual creation idempotent.
 function add_constraint_dual!(
     container::OptimizationContainer,
     sys::IS.InfrastructureSystemsContainer,
     model::ServiceModel{T, D},
 ) where {T <: IS.InfrastructureSystemsComponent, D <: AbstractServiceFormulation}
     if !isempty(get_duals(model))
-        service = get_available_components(model, sys)
+        time_steps = get_time_steps(container)
         for constraint_type in get_duals(model)
-            assign_dual_variable!(container, constraint_type, service, D)
+            for key in _existing_constraint_keys(container, constraint_type, T)
+                # `_existing_constraint_keys` filters on `get_entry_type(key) === constraint_type`,
+                # so `constraint_type` is exactly the key's entry type (avoids re-deriving it).
+                dual_key = ConstraintKey(constraint_type, T, key.meta)
+                haskey(get_duals(container), dual_key) && continue
+                existing = get_constraint(container, key)
+                _assign_dual_from_existing!(container, key, existing, T, time_steps)
+            end
         end
     end
     return
 end
 
-# service formulation
-function assign_dual_variable!(
-    container::OptimizationContainer,
-    constraint_type::Type{<:ConstraintType},
-    service::D,
-    ::Type{<:AbstractServiceFormulation},
-) where {D <: IS.InfrastructureSystemsComponent}
-    time_steps = get_time_steps(container)
-    service_name = IS.get_name(service)
-    add_dual_container!(
-        container,
-        constraint_type,
-        D,
-        [service_name],
-        time_steps;
-        meta = service_name,
-    )
-    return
-end
-
 _existing_constraint_keys(
     container::OptimizationContainer,
     ::Type{T},

src/common_models/add_variable.jl

@@ -77,7 +77,14 @@ function add_variables!(
 end
 
 """
-Add variables to the OptimizationContainer for a service.
+Add variables to the OptimizationContainer for a single service and its contributing
+devices.
+
+All services of a given `(VariableType, ServiceType)` share a single sparse container
+keyed by `(service_name, device_name, time)`, rather than one dense container per
+service disambiguated by a `meta = service_name` field. The container is created lazily
+on the first service of a type, and each subsequent call appends that service's slice, so
+separate formulation groups sharing a service type append to the same container.
 """
 function add_service_variables!(
     container::OptimizationContainer,
@@ -94,37 +101,29 @@ function add_service_variables!(
     @assert !isempty(contributing_devices)
     time_steps = get_time_steps(container)
     settings = get_settings(container)
-
     binary = get_variable_binary(T, U, F)
-
-    variable = add_variable_container!(
-        container,
-        T,
-        U,
-        IS.get_name(service),
-        [IS.get_name(d) for d in contributing_devices],
-        time_steps,
+    s_name = IS.get_name(service)
+    device_names = [IS.get_name(d) for d in contributing_devices]
+    variable = lazy_container_addition!(
+        container, T, U, [s_name], device_names, time_steps; sparse = true,
     )
-
+    jump_model = get_jump_model(container)
     for t in time_steps, d in contributing_devices
         name = IS.get_name(d)
-        variable[name, t] = JuMP.@variable(
-            get_jump_model(container),
-            base_name = "$(T)_$(U)_$(IS.get_name(service))_{$(name), $(t)}",
-            binary = binary
+        var = JuMP.@variable(
+            jump_model,
+            base_name = "$(T)_$(U)_{$(s_name), $(name), $(t)}",
+            binary = binary,
         )
-
+        variable[(s_name, name, t)] = var
         ub = get_variable_upper_bound(T, service, d, F)
-        ub !== nothing && JuMP.set_upper_bound(variable[name, t], ub)
-
+        ub !== nothing && JuMP.set_upper_bound(var, ub)
         lb = get_variable_lower_bound(T, service, d, F)
-        lb !== nothing && !binary && JuMP.set_lower_bound(variable[name, t], lb)
-
+        lb !== nothing && !binary && JuMP.set_lower_bound(var, lb)
         if get_warm_start(settings)
             init = get_variable_warm_start_value(T, d, F)
-            init !== nothing && JuMP.set_start_value(variable[name, t], init)
+            init !== nothing && JuMP.set_start_value(var, init)
         end
     end
-
     return
 end

src/core/definitions.jl

@@ -127,7 +127,6 @@ const M_VALUE = 1e6
 # λ = 0.5 is optimal per Beach et al. (2024), Remark 1.
 const DNMDT_LAMBDA = 0.5
 
-const NO_SERVICE_NAME_PROVIDED = ""
 const UPPER_BOUND = "ub"
 const LOWER_BOUND = "lb"
 const MAX_OPTIMIZE_TRIES = 2

src/core/optimization_container.jl

@@ -731,9 +731,18 @@ function add_variable_container!(
     return _add_container!(container, T, U, JuMP.VariableRef, sparse, axs...; meta = meta)
 end
 
-function _get_pwl_variables_container()
-    contents = Dict{Tuple{String, Int, Int}, JuMP.VariableRef}()
-    return SparseAxisArray(contents)
+"""
+Key tuple type for the empty `SparseAxisArray` auto-created for a `SparseVariableType`.
+
+Defaults to the 3D device-offer PWL shape `(device_name, segment, time)`. A variable type
+that needs an extra axis - e.g. a per-service reserve offer keyed
+`(service_name, device_name, segment, time)` - overrides this method to widen the key. Downstream
+packages extend it for their own sparse variable types.
+"""
+sparse_variable_key_type(::Type{<:SparseVariableType}) = Tuple{String, Int, Int}
+
+function _get_pwl_variables_container(::Type{T}) where {T <: SparseVariableType}
+    return SparseAxisArray(Dict{sparse_variable_key_type(T), JuMP.VariableRef}())
 end
 
 function add_variable_container!(
@@ -746,7 +755,7 @@ function add_variable_container!(
     U <: Union{IS.InfrastructureSystemsComponent, IS.InfrastructureSystemsContainer},
 }
     var_key = VariableKey(T, U, meta)
-    _assign_container!(container.variables, var_key, _get_pwl_variables_container())
+    _assign_container!(container.variables, var_key, _get_pwl_variables_container(T))
     return container.variables[var_key]
 end
 

src/core/service_model.jl

@@ -13,10 +13,9 @@ function _check_service_formulation(::Type{D}) where {D}
 end
 
 """
-Establishes the model for a particular service specified by type. The optional
-`service_name` positional argument assigns the model to a service with that name in the
-template. Uses the keyword argument `feedforwards` to enable passing values between
-operation models at simulation time.
+Establishes the model for all services of a particular type. A `ServiceModel` represents
+every service of its type in the system. Uses the keyword argument `feedforwards` to
+enable passing values between operation models at simulation time.
 
 # Arguments
 
@@ -34,30 +33,29 @@ reserves = ServiceModel(PSY.VariableReserve{PSY.ReserveUp}, RangeReserve)
 mutable struct ServiceModel{D <: IS.InfrastructureSystemsComponent, B}
     # Heterogeneous by design: concrete Vector of the abstract type, not a UnionAll field.
     feedforwards::Vector{AbstractAffectFeedforward}
-    service_name::String
     use_slacks::Bool
     duals::Vector{DataType}
     time_series_names::Dict{Type{<:TimeSeriesParameter}, String}
     attributes::Dict{String, Any}
+    # Per service: service name -> device type -> contributing devices.
     contributing_devices_map::Dict{
-        Type{<:IS.InfrastructureSystemsComponent},
-        Vector{<:IS.InfrastructureSystemsComponent},
+        String,
+        Dict{DataType, Vector{<:IS.InfrastructureSystemsComponent}},
     }
     subsystem::Union{Nothing, String}
     # Maps outage UUIDs to monitored components grouped by device type. PNM indexes DF matrices with UUIDs.
     outages::Dict{Base.UUID, Dict{DataType, Set{String}}}
     function ServiceModel(
         ::Type{D},
-        ::Type{B},
-        service_name::String;
+        ::Type{B};
         use_slacks = false,
         feedforwards = Vector{AbstractAffectFeedforward}(),
         duals = Vector{DataType}(),
         time_series_names = get_default_time_series_names(D, B),
         attributes = Dict{String, Any}(),
         contributing_devices_map = Dict{
-            Type{<:IS.InfrastructureSystemsComponent},
-            Vector{<:IS.InfrastructureSystemsComponent},
+            String,
+            Dict{DataType, Vector{<:IS.InfrastructureSystemsComponent}},
         }(),
     ) where {D <: IS.InfrastructureSystemsComponent, B}
         attributes_for_model = get_default_attributes(D, B)
@@ -69,7 +67,6 @@ mutable struct ServiceModel{D <: IS.InfrastructureSystemsComponent, B}
         _check_service_formulation(B)
         new{D, B}(
             convert(Vector{AbstractAffectFeedforward}, feedforwards),
-            service_name,
             use_slacks,
             duals,
             time_series_names,
@@ -88,53 +85,37 @@ get_formulation(
     ::ServiceModel{D, B},
 ) where {D <: IS.InfrastructureSystemsComponent, B} = B
 get_feedforwards(m::ServiceModel) = m.feedforwards
-get_service_name(m::ServiceModel) = m.service_name
 get_use_slacks(m::ServiceModel) = m.use_slacks
 get_duals(m::ServiceModel) = m.duals
 get_time_series_names(m::ServiceModel) = m.time_series_names
 get_attributes(m::ServiceModel) = m.attributes
 get_attribute(m::ServiceModel, key::String) = get(m.attributes, key, nothing)
+# Whole nested map: service name -> device type -> contributing devices.
 get_contributing_devices_map(m::ServiceModel) = m.contributing_devices_map
-get_contributing_devices_map(m::ServiceModel, key) =
-    get(m.contributing_devices_map, key, nothing)
+# Returned for a service with no entry in the map. Callers treat it as read-only (shared).
+const _EMPTY_CONTRIBUTING_DEVICES_MAP =
+    Dict{DataType, Vector{<:IS.InfrastructureSystemsComponent}}()
+# One service's inner `device type -> devices` map (the empty const if the service is absent).
+get_contributing_devices_map(m::ServiceModel, service_name::AbstractString) =
+    get(m.contributing_devices_map, service_name, _EMPTY_CONTRIBUTING_DEVICES_MAP)
+# All contributing devices across ALL services (flatten the nested map).
+# TODO(services stability): flattening across device types yields a Vector whose element
+# type widens to the abstract common ancestor when a service has more than one contributing
+# device type, so downstream builders lose type stability. Revisit by iterating the
+# per-(device type) map groups (each concretely typed) instead of flattening.
 get_contributing_devices(m::ServiceModel) =
-    [z for x in values(m.contributing_devices_map) for z in x]
+    [z for inner in values(m.contributing_devices_map) for x in values(inner) for z in x]
+# One service's contributing devices (flattened Vector).
+# TODO(services stability): same multi-device-type widening as the all-services flatten
+# above; revisit to iterate the concretely-typed per-device-type map groups.
+get_contributing_devices(m::ServiceModel, service_name::AbstractString) =
+    [z for x in values(get_contributing_devices_map(m, service_name)) for z in x]
 get_subsystem(m::ServiceModel) = m.subsystem
 get_outages(m::ServiceModel) = m.outages
 
 set_subsystem!(m::ServiceModel, id::String) = m.subsystem = id
 
-function ServiceModel(
-    service_type::Type{D},
-    formulation_type::Type{B};
-    use_slacks = false,
-    feedforwards = Vector{AbstractAffectFeedforward}(),
-    duals = Vector{DataType}(),
-    time_series_names = get_default_time_series_names(D, B),
-    attributes = get_default_attributes(D, B),
-) where {D <: IS.InfrastructureSystemsComponent, B}
-    # If more attributes are used later, move free form string to const and organize
-    # attributes
-    attributes_for_model = get_default_attributes(D, B)
-    for (k, v) in attributes
-        attributes_for_model[k] = v
-    end
-    if !haskey(attributes_for_model, "aggregated_service_model")
-        push!(attributes_for_model, "aggregated_service_model" => true)
-    end
-    return ServiceModel(
-        service_type,
-        formulation_type,
-        NO_SERVICE_NAME_PROVIDED;
-        use_slacks,
-        feedforwards,
-        duals,
-        time_series_names,
-        attributes = attributes_for_model,
-    )
-end
-
-function set_model!(dict::Dict, key::Tuple{String, Symbol}, model::ServiceModel)
+function set_model!(dict::Dict, key::Symbol, model::ServiceModel)
     if haskey(dict, key)
         @warn "Overwriting $(key) existing model"
     end
@@ -146,6 +127,6 @@ function set_model!(
     dict::Dict,
     model::ServiceModel{D, B},
 ) where {D <: IS.InfrastructureSystemsComponent, B}
-    set_model!(dict, (get_service_name(model), Symbol(D)), model)
+    set_model!(dict, Symbol(D), model)
     return
 end

src/operation/emulation_model_store.jl

@@ -109,6 +109,28 @@ function write_output!(
     return
 end
 
+# Sparse containers (e.g. reserve variables keyed `(service, device, t)`) are stored as
+# 2D dense with the non-time tuple flattened into encoded `"a__b"` columns, matching the
+# storage `initialize_storage!` allocated via `get_column_names_from_axis_array`. Columns
+# are ordered by their encoded string so the flattened matrix aligns with the
+# pre-allocated dataset's rows (`set_value!` copies positionally).
+function write_output!(
+    store::EmulationModelStore,
+    name::Symbol,
+    key::OptimizationContainerKey,
+    index::EmulationModelIndexType,
+    update_timestamp::Dates.DateTime,
+    array::SparseAxisArray{T, N, K},
+) where {T, N, K <: NTuple{N, Any}}
+    tuple_columns = unique!([k[1:(N - 1)] for k in keys(array.data)])
+    sort!(tuple_columns; by = encode_tuple_to_column)
+    matrix = _to_matrix(array, tuple_columns)
+    columns = encode_tuple_to_column.(tuple_columns)
+    dense = DenseAxisArray(permutedims(matrix), columns, 1:size(matrix, 1))
+    write_output!(store, name, key, index, update_timestamp, dense)
+    return
+end
+
 function read_outputs(
     store::EmulationModelStore{InMemoryDataset},
     key::OptimizationContainerKey;

src/operation/problem_template.jl

@@ -1,6 +1,6 @@
 
 const DevicesModelContainer = Dict{Symbol, DeviceModel}
-const ServicesModelContainer = Dict{Tuple{String, Symbol}, ServiceModel}
+const ServicesModelContainer = Dict{Symbol, ServiceModel}
 
 abstract type AbstractProblemTemplate end
 

test/test_emulation_model_store.jl

@@ -23,3 +23,30 @@ DatasetContainer fields (Task 2.10).
     @test isempty(store)
     @test isempty(IOM.list_keys(store, IOM.VariableType))
 end
+
+@testset "EmulationModelStore sparse 3D container write/read round trip" begin
+    # Reserve-style containers are sparse and keyed `(service, device, time)`; the store
+    # flattens the leading dims to encoded `"service__device"` columns. This mirrors the
+    # DecisionModelStore path but the EmulationModelStore write must line the flattened
+    # rows up with the pre-allocated dataset (`set_value!` copies positionally).
+    store = IOM.EmulationModelStore()
+    key = IOM.VariableKey(TestVariableType, MockComponentType)
+    values = Dict(
+        ("s1", "d1", 1) => 1.0,
+        ("s1", "d2", 1) => 2.0,
+        ("s2", "d1", 1) => 3.0,
+    )
+    sparse = JuMP.Containers.SparseAxisArray(values)
+    # Storage is pre-allocated exactly as initialize_storage! would, from the encoded
+    # column names.
+    cols = IOM.get_column_names_from_axis_array(key, sparse)[1]
+    storage = DenseAxisArray(fill(NaN, length(cols), 1), cols, 1:1)
+    IOM.set_dataset!(store.data_container, key, IOM.InMemoryDataset(storage))
+
+    IOM.write_output!(store, :variables, key, 1, Dates.DateTime(2024, 1, 1), sparse)
+    out = IOM.read_outputs(store, key)
+
+    @test out["s1__d1", 1] == 1.0
+    @test out["s1__d2", 1] == 2.0
+    @test out["s2__d1", 1] == 3.0
+end