diff --git a/src/InfrastructureSystems.jl b/src/InfrastructureSystems.jl
index bbf1fa251..b61c5913d 100644
--- a/src/InfrastructureSystems.jl
+++ b/src/InfrastructureSystems.jl
@@ -8,6 +8,14 @@ export PiecewisePointCurve, PiecewiseIncrementalCurve, PiecewiseAverageCurve
 export TimeSeriesLinearCurve, TimeSeriesQuadraticCurve, TimeSeriesPiecewisePointCurve
 export TimeSeriesPiecewiseIncrementalCurve, TimeSeriesPiecewiseAverageCurve
 
+# Units interface: declared here, methods implemented by domain packages
+# (e.g., PowerSystems.jl provides power-domain `get_value`/`set_value` methods).
+"Get a field value with optional unit conversion. Methods are provided by domain packages."
+function get_value end
+"Set a field value with optional unit conversion. Methods are provided by domain packages."
+function set_value end
+export get_value, set_value
+
 import Base: @kwdef
 import CSV
 import DataFrames
@@ -134,6 +142,7 @@ end
 get_internal(value::InfrastructureSystemsComponent) = value.internal
 
 include("common.jl")
+include("relative_units.jl")
 include("random_seed.jl")
 include("utils/timers.jl")
 include("utils/assert_op.jl")
@@ -208,4 +217,5 @@ include("function_data/make_convex.jl")
 include("deprecated.jl")
 include("Optimization/Optimization.jl")
 include("Simulation/Simulation.jl")
+
 end # module
diff --git a/src/cost_aliases.jl b/src/cost_aliases.jl
index c42144322..4faabf83a 100644
--- a/src/cost_aliases.jl
+++ b/src/cost_aliases.jl
@@ -8,8 +8,8 @@
 # methods being defined for all the `ValueCurve{FunctionData}` types, not just the ones we
 # have here nicely packaged and presented to the user.
 
-# Default `is_cost_alias` is defined in value_curve.jl so it's available to
-# time_series_value_curve.jl show methods (included before this file).
+"Whether there is a cost alias for the instance or type under consideration"
+is_cost_alias(::Union{ValueCurve, Type{<:ValueCurve}}) = false
 
 """
     LinearCurve(proportional_term::Float64)
@@ -29,6 +29,7 @@ curve = LinearCurve(50.0, 100.0)
 const LinearCurve = InputOutputCurve{LinearFunctionData}
 
 is_cost_alias(::Union{LinearCurve, Type{LinearCurve}}) = true
+simple_type_name(::LinearCurve) = "LinearCurve"
 
 InputOutputCurve{LinearFunctionData}(proportional_term::Real) =
     InputOutputCurve(LinearFunctionData(proportional_term))
@@ -44,7 +45,7 @@ get_constant_term(vc::LinearCurve) = get_constant_term(get_function_data(vc))
 
 Base.show(io::IO, vc::LinearCurve) =
     if isnothing(get_input_at_zero(vc))
-        print(io, "$(typeof(vc))($(get_proportional_term(vc)), $(get_constant_term(vc)))")
+        print(io, "LinearCurve($(get_proportional_term(vc)), $(get_constant_term(vc)))")
     else
         Base.show_default(io, vc)
     end
@@ -67,6 +68,7 @@ curve = QuadraticCurve(0.002, 25.0, 150.0)
 const QuadraticCurve = InputOutputCurve{QuadraticFunctionData}
 
 is_cost_alias(::Union{QuadraticCurve, Type{QuadraticCurve}}) = true
+simple_type_name(::QuadraticCurve) = "QuadraticCurve"
 
 InputOutputCurve{QuadraticFunctionData}(quadratic_term, proportional_term, constant_term) =
     InputOutputCurve(
@@ -86,7 +88,7 @@ Base.show(io::IO, vc::QuadraticCurve) =
     if isnothing(get_input_at_zero(vc))
         print(
             io,
-            "$(typeof(vc))($(get_quadratic_term(vc)), $(get_proportional_term(vc)), $(get_constant_term(vc)))",
+            "QuadraticCurve($(get_quadratic_term(vc)), $(get_proportional_term(vc)), $(get_constant_term(vc)))",
         )
     else
         Base.show_default(io, vc)
@@ -112,6 +114,7 @@ curve = PiecewisePointCurve([(100.0, 400.0), (200.0, 900.0), (300.0, 1500.0)])
 const PiecewisePointCurve = InputOutputCurve{PiecewiseLinearData}
 
 is_cost_alias(::Union{PiecewisePointCurve, Type{PiecewisePointCurve}}) = true
+simple_type_name(::PiecewisePointCurve) = "PiecewisePointCurve"
 
 InputOutputCurve{PiecewiseLinearData}(points::Vector) =
     InputOutputCurve(PiecewiseLinearData(points))
@@ -131,7 +134,7 @@ get_slopes(vc::PiecewisePointCurve) = get_slopes(get_function_data(vc))
 # Here we manually circumvent the @NamedTuple{x::Float64, y::Float64} type annotation, but we keep things looking like named tuples
 Base.show(io::IO, vc::PiecewisePointCurve) =
     if isnothing(get_input_at_zero(vc))
-        print(io, "$(typeof(vc))([$(join(get_points(vc), ", "))])")
+        print(io, "PiecewisePointCurve([$(join(get_points(vc), ", "))])")
     else
         Base.show_default(io, vc)
     end
@@ -159,6 +162,7 @@ curve = PiecewiseIncrementalCurve(500.0, [100.0, 150.0, 200.0], [30.0, 35.0])
 const PiecewiseIncrementalCurve = IncrementalCurve{PiecewiseStepData}
 
 is_cost_alias(::Union{PiecewiseIncrementalCurve, Type{PiecewiseIncrementalCurve}}) = true
+simple_type_name(::PiecewiseIncrementalCurve) = "PiecewiseIncrementalCurve"
 
 IncrementalCurve{PiecewiseStepData}(initial_input, x_coords::Vector, slopes::Vector) =
     IncrementalCurve(PiecewiseStepData(x_coords, slopes), initial_input)
@@ -181,9 +185,9 @@ Base.show(io::IO, vc::PiecewiseIncrementalCurve) =
     print(
         io,
         if isnothing(get_input_at_zero(vc))
-            "$(typeof(vc))($(get_initial_input(vc)), $(get_x_coords(vc)), $(get_slopes(vc)))"
+            "PiecewiseIncrementalCurve($(get_initial_input(vc)), $(get_x_coords(vc)), $(get_slopes(vc)))"
         else
-            "$(typeof(vc))($(get_input_at_zero(vc)), $(get_initial_input(vc)), $(get_x_coords(vc)), $(get_slopes(vc)))"
+            "PiecewiseIncrementalCurve($(get_input_at_zero(vc)), $(get_initial_input(vc)), $(get_x_coords(vc)), $(get_slopes(vc)))"
         end,
     )
 
@@ -202,6 +206,7 @@ input). If your data gives incremental/marginal rates instead, use
 const PiecewiseAverageCurve = AverageRateCurve{PiecewiseStepData}
 
 is_cost_alias(::Union{PiecewiseAverageCurve, Type{PiecewiseAverageCurve}}) = true
+simple_type_name(::PiecewiseAverageCurve) = "PiecewiseAverageCurve"
 
 AverageRateCurve{PiecewiseStepData}(initial_input, x_coords::Vector, y_coords::Vector) =
     AverageRateCurve(PiecewiseStepData(x_coords, y_coords), initial_input)
@@ -216,7 +221,7 @@ Base.show(io::IO, vc::PiecewiseAverageCurve) =
     if isnothing(get_input_at_zero(vc))
         print(
             io,
-            "$(typeof(vc))($(get_initial_input(vc)), $(get_x_coords(vc)), $(get_average_rates(vc)))",
+            "PiecewiseAverageCurve($(get_initial_input(vc)), $(get_x_coords(vc)), $(get_average_rates(vc)))",
         )
     else
         Base.show_default(io, vc)
@@ -237,6 +242,7 @@ const TimeSeriesLinearCurve =
     TimeSeriesInputOutputCurve{TimeSeriesFunctionData{LinearFunctionData}}
 
 is_cost_alias(::Union{TimeSeriesLinearCurve, Type{TimeSeriesLinearCurve}}) = true
+simple_type_name(::TimeSeriesLinearCurve) = "TimeSeriesLinearCurve"
 
 TimeSeriesInputOutputCurve{TimeSeriesFunctionData{LinearFunctionData}}(
     key::TimeSeriesKey,
@@ -244,7 +250,7 @@ TimeSeriesInputOutputCurve{TimeSeriesFunctionData{LinearFunctionData}}(
 
 Base.show(io::IO, vc::TimeSeriesLinearCurve) =
     if isnothing(get_input_at_zero(vc))
-        print(io, "$(typeof(vc))($(_ts_key_repr(get_time_series_key(vc))))")
+        print(io, "$(simple_type_name(vc))($(_ts_key_repr(get_time_series_key(vc))))")
     else
         Base.show_default(io, vc)
     end
@@ -259,6 +265,7 @@ const TimeSeriesQuadraticCurve =
     TimeSeriesInputOutputCurve{TimeSeriesFunctionData{QuadraticFunctionData}}
 
 is_cost_alias(::Union{TimeSeriesQuadraticCurve, Type{TimeSeriesQuadraticCurve}}) = true
+simple_type_name(::TimeSeriesQuadraticCurve) = "TimeSeriesQuadraticCurve"
 
 TimeSeriesInputOutputCurve{TimeSeriesFunctionData{QuadraticFunctionData}}(
     key::TimeSeriesKey,
@@ -266,7 +273,7 @@ TimeSeriesInputOutputCurve{TimeSeriesFunctionData{QuadraticFunctionData}}(
 
 Base.show(io::IO, vc::TimeSeriesQuadraticCurve) =
     if isnothing(get_input_at_zero(vc))
-        print(io, "$(typeof(vc))($(_ts_key_repr(get_time_series_key(vc))))")
+        print(io, "$(simple_type_name(vc))($(_ts_key_repr(get_time_series_key(vc))))")
     else
         Base.show_default(io, vc)
     end
@@ -283,6 +290,7 @@ const TimeSeriesPiecewisePointCurve =
 is_cost_alias(
     ::Union{TimeSeriesPiecewisePointCurve, Type{TimeSeriesPiecewisePointCurve}},
 ) = true
+simple_type_name(::TimeSeriesPiecewisePointCurve) = "TimeSeriesPiecewisePointCurve"
 
 TimeSeriesInputOutputCurve{TimeSeriesFunctionData{PiecewiseLinearData}}(
     key::TimeSeriesKey,
@@ -290,7 +298,7 @@ TimeSeriesInputOutputCurve{TimeSeriesFunctionData{PiecewiseLinearData}}(
 
 Base.show(io::IO, vc::TimeSeriesPiecewisePointCurve) =
     if isnothing(get_input_at_zero(vc))
-        print(io, "$(typeof(vc))($(_ts_key_repr(get_time_series_key(vc))))")
+        print(io, "$(simple_type_name(vc))($(_ts_key_repr(get_time_series_key(vc))))")
     else
         Base.show_default(io, vc)
     end
@@ -310,6 +318,8 @@ is_cost_alias(
         Type{TimeSeriesPiecewiseIncrementalCurve},
     },
 ) = true
+simple_type_name(::TimeSeriesPiecewiseIncrementalCurve) =
+    "TimeSeriesPiecewiseIncrementalCurve"
 
 TimeSeriesIncrementalCurve{TimeSeriesFunctionData{PiecewiseStepData}}(
     key::TimeSeriesKey,
@@ -330,9 +340,9 @@ Base.show(io::IO, vc::TimeSeriesPiecewiseIncrementalCurve) =
     print(
         io,
         if isnothing(get_input_at_zero(vc))
-            "$(typeof(vc))($(_ts_key_repr(get_time_series_key(vc))), $(_ts_key_repr(get_initial_input(vc))))"
+            "$(simple_type_name(vc))($(_ts_key_repr(get_time_series_key(vc))), $(_ts_key_repr(get_initial_input(vc))))"
         else
-            "$(typeof(vc))($(_ts_key_repr(get_time_series_key(vc))), $(_ts_key_repr(get_initial_input(vc))), $(_ts_key_repr(get_input_at_zero(vc))))"
+            "$(simple_type_name(vc))($(_ts_key_repr(get_time_series_key(vc))), $(_ts_key_repr(get_initial_input(vc))), $(_ts_key_repr(get_input_at_zero(vc))))"
         end,
     )
 
@@ -351,6 +361,7 @@ is_cost_alias(
         Type{TimeSeriesPiecewiseAverageCurve},
     },
 ) = true
+simple_type_name(::TimeSeriesPiecewiseAverageCurve) = "TimeSeriesPiecewiseAverageCurve"
 
 TimeSeriesAverageRateCurve{TimeSeriesFunctionData{PiecewiseStepData}}(
     key::TimeSeriesKey,
@@ -371,7 +382,7 @@ Base.show(io::IO, vc::TimeSeriesPiecewiseAverageCurve) =
     if isnothing(get_input_at_zero(vc))
         print(
             io,
-            "$(typeof(vc))($(_ts_key_repr(get_time_series_key(vc))), $(_ts_key_repr(get_initial_input(vc))))",
+            "$(simple_type_name(vc))($(_ts_key_repr(get_time_series_key(vc))), $(_ts_key_repr(get_initial_input(vc))))",
         )
     else
         Base.show_default(io, vc)
diff --git a/src/production_variable_cost_curve.jl b/src/production_variable_cost_curve.jl
index 154457ce1..89e21234d 100644
--- a/src/production_variable_cost_curve.jl
+++ b/src/production_variable_cost_curve.jl
@@ -1,21 +1,34 @@
 """
-Supertype for production variable cost curve representations, parameterized by
-a [`ValueCurve`](@ref) type.
+Supertype for production variable cost curve representations.
+
+Parameterized by a [`ValueCurve`](@ref) type `T` and an
+[`AbstractUnitSystem`](@ref) type `U`. `U` is a compile-time marker for the
+`power_units` of the x-axis; it replaces the old `power_units::UnitSystem`
+runtime field. This lets unit-dependent operations dispatch directly on the
+type parameter rather than going through a stateful runtime check plus
+`Val`-wrapping, eliminating a class of type instabilities downstream.
 
 Concrete subtypes include [`CostCurve`](@ref) and [`FuelCurve`](@ref).
 """
-abstract type ProductionVariableCostCurve{T <: ValueCurve} end
+abstract type ProductionVariableCostCurve{T <: ValueCurve, U <: AbstractUnitSystem} end
 
-serialize(val::ProductionVariableCostCurve) = serialize_struct(val)
-deserialize(T::Type{<:ProductionVariableCostCurve}, val::Dict) =
-    deserialize_struct(T, val)
+"""
+`ProductionVariableCostCurve{T}` with any unit system. Use at `isa` sites or in
+docstrings where the curve's unit-system parameter doesn't matter.
+"""
+const AnyProductionVariableCostCurve{T} =
+    ProductionVariableCostCurve{T, U} where {U <: AbstractUnitSystem}
 
 "Get the underlying `ValueCurve` representation of this `ProductionVariableCostCurve`"
 get_value_curve(cost::ProductionVariableCostCurve) = cost.value_curve
 "Get the variable operation and maintenance cost in currency/(power_units h)"
 get_vom_cost(cost::ProductionVariableCostCurve) = cost.vom_cost
-"Get the units for the x-axis of the curve"
-get_power_units(cost::ProductionVariableCostCurve) = cost.power_units
+"""
+Get the units marker for the x-axis of the curve as an instance of the
+second type parameter (e.g. `NaturalUnit()`, `SystemBaseUnit()`,
+`DeviceBaseUnit()`).
+"""
+get_power_units(::ProductionVariableCostCurve{T, U}) where {T, U} = U()
 "Get the `FunctionData` representation of this `ProductionVariableCostCurve`'s `ValueCurve`"
 get_function_data(cost::ProductionVariableCostCurve) =
     get_function_data(get_value_curve(cost))
@@ -57,120 +70,181 @@ Base.hash(a::ProductionVariableCostCurve, h::UInt) = hash_from_fields(a, h)
 $(TYPEDEF)
 $(TYPEDFIELDS)
 
+    CostCurve(value_curve)
+    CostCurve(value_curve, power_units)
+    CostCurve(value_curve, vom_cost)
     CostCurve(value_curve, power_units, vom_cost)
     CostCurve(; value_curve, power_units, vom_cost)
 
 Direct representation of the variable operation cost of a power plant in currency. Composed
 of a [`ValueCurve`](@ref) that may represent input-output, incremental, or average rate
-data. The default units for the x-axis are MW and can be specified with
-`power_units`.
+data. The x-axis units are encoded as the second type parameter `U <: AbstractUnitSystem`;
+`power_units` at construction is the singleton instance `U()` (default `NaturalUnit()`).
 """
-@kwdef struct CostCurve{T <: ValueCurve} <: ProductionVariableCostCurve{T}
+struct CostCurve{T <: ValueCurve, U <: AbstractUnitSystem} <:
+       ProductionVariableCostCurve{T, U}
     "The underlying `ValueCurve` representation of this `ProductionVariableCostCurve`"
     value_curve::T
-    "(default: natural units (MW)) The units for the x-axis of the curve"
-    power_units::UnitSystem = UnitSystem.NATURAL_UNITS
     "(default of 0) Additional proportional Variable Operation and Maintenance Cost in
     \$/(power_unit h), represented as a [`LinearCurve`](@ref)"
-    vom_cost::LinearCurve = LinearCurve(0.0)
+    vom_cost::LinearCurve
+
+    CostCurve{T, U}(value_curve::T, vom_cost::LinearCurve) where {T, U} =
+        new{T, U}(value_curve, vom_cost)
 end
 
-CostCurve(value_curve) = CostCurve(; value_curve)
-CostCurve(value_curve, vom_cost::LinearCurve) =
-    CostCurve(; value_curve, vom_cost = vom_cost)
-CostCurve(value_curve, power_units::UnitSystem) =
-    CostCurve(; value_curve, power_units = power_units)
+CostCurve{T, U}(;
+    value_curve::T,
+    vom_cost::LinearCurve = LinearCurve(0.0),
+) where {T, U} = CostCurve{T, U}(value_curve, vom_cost)
+
+# Outer constructors — default U = NaturalUnit when not specified
+CostCurve(value_curve::T) where {T <: ValueCurve} =
+    CostCurve{T, NaturalUnit}(; value_curve)
+CostCurve(value_curve::T, vom_cost::LinearCurve) where {T <: ValueCurve} =
+    CostCurve{T, NaturalUnit}(; value_curve, vom_cost)
+CostCurve(
+    value_curve::T,
+    power_units::U,
+) where {T <: ValueCurve, U <: AbstractUnitSystem} =
+    CostCurve{T, U}(; value_curve)
+CostCurve(
+    value_curve::T,
+    power_units::U,
+    vom_cost::LinearCurve,
+) where {T <: ValueCurve, U <: AbstractUnitSystem} =
+    CostCurve{T, U}(; value_curve, vom_cost)
 
-Base.:(==)(a::CostCurve, b::CostCurve) =
-    (get_value_curve(a) == get_value_curve(b)) &&
-    (get_power_units(a) == get_power_units(b)) &&
-    (get_vom_cost(a) == get_vom_cost(b))
+# Keyword-based constructor exposing `power_units`, replacing the former field default
+function CostCurve(;
+    value_curve,
+    power_units::AbstractUnitSystem = NaturalUnit(),
+    vom_cost::LinearCurve = LinearCurve(0.0),
+)
+    return CostCurve{typeof(value_curve), typeof(power_units)}(; value_curve, vom_cost)
+end
 
 "Get a `CostCurve` representing zero variable cost"
 Base.zero(::Union{CostCurve, Type{CostCurve}}) = CostCurve(zero(ValueCurve))
 
+"""
+`CostCurve{T}` with any unit system. Equivalent to `CostCurve{T, U} where U`;
+use at `isa` sites where the unit-system parameter doesn't matter.
+"""
+const AnyCostCurve{T} = CostCurve{T, U} where {U <: AbstractUnitSystem}
+
 """
 $(TYPEDEF)
 $(TYPEDFIELDS)
 
-    FuelCurve(value_curve, power_units, fuel_cost, startup_fuel_offtake, vom_cost)
     FuelCurve(value_curve, fuel_cost)
     FuelCurve(value_curve, fuel_cost, startup_fuel_offtake, vom_cost)
     FuelCurve(value_curve, power_units, fuel_cost)
+    FuelCurve(value_curve, power_units, fuel_cost, startup_fuel_offtake, vom_cost)
     FuelCurve(; value_curve, power_units, fuel_cost, startup_fuel_offtake, vom_cost)
 
 Representation of the variable operation cost of a power plant in terms of fuel (MBTU,
 liters, m^3, etc.), coupled with a conversion factor between fuel and currency. Composed of
 a [`ValueCurve`](@ref) that may represent input-output, incremental, or average rate data.
-The default units for the x-axis are MW and can be specified with `power_units`.
+The x-axis units are encoded as the second type parameter `U <: AbstractUnitSystem`;
+`power_units` at construction is the singleton instance `U()` (default `NaturalUnit()`).
 """
-@kwdef struct FuelCurve{T <: ValueCurve} <: ProductionVariableCostCurve{T}
+struct FuelCurve{T <: ValueCurve, U <: AbstractUnitSystem} <:
+       ProductionVariableCostCurve{T, U}
     "The underlying `ValueCurve` representation of this `ProductionVariableCostCurve`"
     value_curve::T
-    "(default: natural units (MW)) The units for the x-axis of the curve"
-    power_units::UnitSystem = UnitSystem.NATURAL_UNITS
     "Either a fixed value for fuel cost or the [`TimeSeriesKey`](@ref) to a fuel cost time series"
     fuel_cost::Union{Float64, TimeSeriesKey}
     "(default of 0) Fuel consumption at the unit startup proceedure. Additional cost to the startup costs and related only to the initial fuel required to start the unit.
     represented as a [`LinearCurve`](@ref)"
-    startup_fuel_offtake::LinearCurve = LinearCurve(0.0)
+    startup_fuel_offtake::LinearCurve
     "(default of 0) Additional proportional Variable Operation and Maintenance Cost in \$/(power_unit h)
     represented as a [`LinearCurve`](@ref)"
-    vom_cost::LinearCurve = LinearCurve(0.0)
+    vom_cost::LinearCurve
+
+    FuelCurve{T, U}(
+        value_curve::T,
+        fuel_cost::Union{Float64, TimeSeriesKey},
+        startup_fuel_offtake::LinearCurve,
+        vom_cost::LinearCurve,
+    ) where {T, U} =
+        new{T, U}(value_curve, fuel_cost, startup_fuel_offtake, vom_cost)
 end
 
-function FuelCurve(
-    value_curve::ValueCurve,
-    power_units::UnitSystem,
-    fuel_cost::Real,
+FuelCurve{T, U}(;
+    value_curve::T,
+    fuel_cost::Union{Float64, TimeSeriesKey},
+    startup_fuel_offtake::LinearCurve = LinearCurve(0.0),
+    vom_cost::LinearCurve = LinearCurve(0.0),
+) where {T, U} =
+    FuelCurve{T, U}(value_curve, fuel_cost, startup_fuel_offtake, vom_cost)
+
+# Outer constructors — mirror the CostCurve style
+FuelCurve(value_curve::T, fuel_cost::Real) where {T <: ValueCurve} =
+    FuelCurve{T, NaturalUnit}(; value_curve, fuel_cost = Float64(fuel_cost))
+FuelCurve(value_curve::T, fuel_cost::TimeSeriesKey) where {T <: ValueCurve} =
+    FuelCurve{T, NaturalUnit}(; value_curve, fuel_cost)
+
+FuelCurve(
+    value_curve::T,
+    fuel_cost::Union{Real, TimeSeriesKey},
     startup_fuel_offtake::LinearCurve,
     vom_cost::LinearCurve,
+) where {T <: ValueCurve} = FuelCurve{T, NaturalUnit}(;
+    value_curve,
+    fuel_cost = fuel_cost isa Real ? Float64(fuel_cost) : fuel_cost,
+    startup_fuel_offtake,
+    vom_cost,
 )
-    return FuelCurve(
-        value_curve,
-        power_units,
-        Float64(fuel_cost),
-        startup_fuel_offtake,
-        vom_cost,
-    )
-end
-
-function FuelCurve(value_curve, fuel_cost)
-    FuelCurve(; value_curve, fuel_cost)
-end
 
-function FuelCurve(
+FuelCurve(
+    value_curve::T,
+    power_units::U,
+    fuel_cost::Union{Real, TimeSeriesKey},
+) where {T <: ValueCurve, U <: AbstractUnitSystem} = FuelCurve{T, U}(;
     value_curve,
-    fuel_cost::Union{Float64, TimeSeriesKey},
+    fuel_cost = fuel_cost isa Real ? Float64(fuel_cost) : fuel_cost,
+)
+
+FuelCurve(
+    value_curve::T,
+    power_units::U,
+    fuel_cost::Union{Real, TimeSeriesKey},
     startup_fuel_offtake::LinearCurve,
     vom_cost::LinearCurve,
+) where {T <: ValueCurve, U <: AbstractUnitSystem} = FuelCurve{T, U}(;
+    value_curve,
+    fuel_cost = fuel_cost isa Real ? Float64(fuel_cost) : fuel_cost,
+    startup_fuel_offtake,
+    vom_cost,
 )
-    return FuelCurve(;
-        value_curve,
-        fuel_cost,
-        startup_fuel_offtake = startup_fuel_offtake,
-        vom_cost = vom_cost,
-    )
-end
 
-function FuelCurve(
+# Keyword-based constructor exposing `power_units`
+function FuelCurve(;
     value_curve,
-    power_units::UnitSystem,
-    fuel_cost::Union{Float64, TimeSeriesKey},
+    power_units::AbstractUnitSystem = NaturalUnit(),
+    fuel_cost::Union{Real, TimeSeriesKey},
+    startup_fuel_offtake::LinearCurve = LinearCurve(0.0),
+    vom_cost::LinearCurve = LinearCurve(0.0),
 )
-    FuelCurve(; value_curve, power_units = power_units, fuel_cost = fuel_cost)
+    fc = fuel_cost isa Real ? Float64(fuel_cost) : fuel_cost
+    return FuelCurve{typeof(value_curve), typeof(power_units)}(;
+        value_curve,
+        fuel_cost = fc,
+        startup_fuel_offtake,
+        vom_cost,
+    )
 end
 
-Base.:(==)(a::FuelCurve, b::FuelCurve) =
-    (get_value_curve(a) == get_value_curve(b)) &&
-    (get_power_units(a) == get_power_units(b)) &&
-    (get_fuel_cost(a) == get_fuel_cost(b)) &&
-    (get_startup_fuel_offtake(a) == get_startup_fuel_offtake(b)) &&
-    (get_vom_cost(a) == get_vom_cost(b))
-
 "Get a `FuelCurve` representing zero fuel usage and zero fuel cost"
 Base.zero(::Union{FuelCurve, Type{FuelCurve}}) = FuelCurve(zero(ValueCurve), 0.0)
 
+"""
+`FuelCurve{T}` with any unit system. Equivalent to `FuelCurve{T, U} where U`;
+use at `isa` sites where the unit-system parameter doesn't matter.
+"""
+const AnyFuelCurve{T} = FuelCurve{T, U} where {U <: AbstractUnitSystem}
+
 "Get the fuel cost or the name of the fuel cost time series"
 get_fuel_cost(cost::FuelCurve) = cost.fuel_cost
 "Get the function for the fuel consumption at startup"
@@ -186,6 +260,59 @@ is_time_series_backed(cost::FuelCurve) =
     is_time_series_backed(get_value_curve(cost)) ||
     is_time_series_backed(get_fuel_cost(cost))
 
+# ── Serialization ─────────────────────────────────────────────────────────────
+# The U type parameter has no corresponding field, so we serialize it under the
+# conventional "power_units" key (preserving the field name from the previous
+# schema) and reconstruct it at deserialize time.
+
+# Accept both new (type-name) and legacy (UnitSystem enum-value-name) encodings so
+# existing serialized fixtures keep deserializing.
+_unit_system_instance(name::AbstractString) = _unit_system_instance(String(name))
+function _unit_system_instance(name::String)
+    if name in ("NATURAL_UNITS", "NaturalUnit")
+        return NaturalUnit()
+    elseif name in ("SYSTEM_BASE", "SystemBaseUnit")
+        return SystemBaseUnit()
+    elseif name in ("DEVICE_BASE", "DeviceBaseUnit")
+        return DeviceBaseUnit()
+    end
+    T = getproperty(@__MODULE__, Symbol(name))
+    T <: AbstractUnitSystem ||
+        throw(ArgumentError("$name is not a subtype of AbstractUnitSystem"))
+    return T()
+end
+
+function serialize(val::ProductionVariableCostCurve)
+    data = serialize_struct(val)
+    data["power_units"] = string(nameof(typeof(get_power_units(val))))
+    return data
+end
+
+function deserialize(::Type{CostCurve}, data::Dict)
+    vc_data = data["value_curve"]
+    vc_type = get_type_from_serialization_data(vc_data)
+    value_curve = deserialize(vc_type, vc_data)
+    vom_cost = deserialize(LinearCurve, data["vom_cost"])
+    power_units = _unit_system_instance(data["power_units"])
+    return CostCurve(value_curve, power_units, vom_cost)
+end
+
+function deserialize(::Type{FuelCurve}, data::Dict)
+    vc_data = data["value_curve"]
+    vc_type = get_type_from_serialization_data(vc_data)
+    value_curve = deserialize(vc_type, vc_data)
+    startup = deserialize(LinearCurve, data["startup_fuel_offtake"])
+    vom = deserialize(LinearCurve, data["vom_cost"])
+    fuel_cost_raw = data["fuel_cost"]
+    fuel_cost = if fuel_cost_raw isa Dict
+        deserialize(TimeSeriesKey, fuel_cost_raw)
+    else
+        Float64(fuel_cost_raw)
+    end
+    power_units = _unit_system_instance(data["power_units"])
+    return FuelCurve(value_curve, power_units, fuel_cost, startup, vom)
+end
+
 Base.show(io::IO, m::MIME"text/plain", curve::ProductionVariableCostCurve) =
     (get(io, :compact, false)::Bool ? _show_compact : _show_expanded)(io, m, curve)
 
@@ -193,7 +320,7 @@ Base.show(io::IO, m::MIME"text/plain", curve::ProductionVariableCostCurve) =
 function _show_compact(io::IO, ::MIME"text/plain", curve::CostCurve)
     print(
         io,
-        "$(nameof(typeof(curve))) with power_units $(curve.power_units), vom_cost $(curve.vom_cost), and value_curve:\n  ",
+        "$(nameof(typeof(curve))) with power_units $(get_power_units(curve)), vom_cost $(curve.vom_cost), and value_curve:\n  ",
     )
     vc_printout = sprint(show, "text/plain", curve.value_curve; context = io)  # Capture the value_curve `show` so we can indent it
     print(io, replace(vc_printout, "\n" => "\n  "))
@@ -202,7 +329,7 @@ end
 function _show_compact(io::IO, ::MIME"text/plain", curve::FuelCurve)
     print(
         io,
-        "$(nameof(typeof(curve))) with power_units $(curve.power_units), fuel_cost $(curve.fuel_cost), startup_fuel_offtake $(curve.startup_fuel_offtake), vom_cost $(curve.vom_cost), and value_curve:\n  ",
+        "$(nameof(typeof(curve))) with power_units $(get_power_units(curve)), fuel_cost $(curve.fuel_cost), startup_fuel_offtake $(curve.startup_fuel_offtake), vom_cost $(curve.vom_cost), and value_curve:\n  ",
     )
     vc_printout = sprint(show, "text/plain", curve.value_curve; context = io)
     print(io, replace(vc_printout, "\n" => "\n  "))
@@ -216,4 +343,6 @@ function _show_expanded(io::IO, ::MIME"text/plain", curve::ProductionVariableCos
             replace(sprint(show, "text/plain", val; context = io), "\n" => "\n  ")
         print(io, "\n  $(field_name): $val_printout")
     end
+    # Surface the type-parameter `power_units` even though it isn't a field
+    print(io, "\n  power_units: $(get_power_units(curve))")
 end
diff --git a/src/relative_units.jl b/src/relative_units.jl
new file mode 100644
index 000000000..9ea59f70e
--- /dev/null
+++ b/src/relative_units.jl
@@ -0,0 +1,165 @@
+###############################
+# Relative (per-unit) markers and RelativeQuantity wrapper.
+#
+# These types are domain-agnostic — they express "device base" / "system base"
+# / "natural unit" without assuming any particular physical domain. Downstream
+# packages (e.g. PowerSystems) attach domain-specific meaning via categories
+# and conversions.
+###############################
+
+"""
+Supertype for all unit-system markers (relative and natural). Used as the
+`U` type parameter on `ProductionVariableCostCurve` and related parametric
+types so that the unit system can be dispatched on at compile time.
+"""
+abstract type AbstractUnitSystem end
+
+"""
+Supertype of per-unit (relative) unit markers.
+"""
+abstract type AbstractRelativeUnit <: AbstractUnitSystem end
+
+"""
+Device base per-unit. Values are normalized to the component's own base.
+"""
+struct DeviceBaseUnit <: AbstractRelativeUnit end
+
+"""
+System base per-unit. Values are normalized to the system's base.
+"""
+struct SystemBaseUnit <: AbstractRelativeUnit end
+
+"""
+Natural units. When used as a target, returns the value with the
+domain-appropriate unit attached (e.g. MW for power, Ω for impedance).
+Deliberately *not* `<: AbstractRelativeUnit` — "convert to NU" yields a
+`Unitful.Quantity`, not a `RelativeQuantity` — but it is a peer under
+`AbstractUnitSystem`.
+"""
+struct NaturalUnit <: AbstractUnitSystem end
+
+const DU = DeviceBaseUnit()
+const SU = SystemBaseUnit()
+const NU = NaturalUnit()
+
+"""
+    RelativeQuantity{T<:Number, U<:AbstractRelativeUnit} <: Number
+
+A quantity tagged with a per-unit marker.
+
+# Examples
+```julia
+0.6 * DU  # 0.6 per-unit on device base
+0.3 * SU  # 0.3 per-unit on system base
+```
+"""
+struct RelativeQuantity{T <: Number, U <: AbstractRelativeUnit} <: Number
+    value::T
+    unit::U
+end
+
+# Construction via multiplication
+Base.:*(a::Number, b::AbstractRelativeUnit) = RelativeQuantity(a, b)
+Base.:*(b::AbstractRelativeUnit, a::Number) = RelativeQuantity(a, b)
+
+# Arithmetic — same unit type only
+Base.:+(a::RelativeQuantity{T, U}, b::RelativeQuantity{S, U}) where {T, S, U} =
+    RelativeQuantity(a.value + b.value, a.unit)
+Base.:-(a::RelativeQuantity{T, U}, b::RelativeQuantity{S, U}) where {T, S, U} =
+    RelativeQuantity(a.value - b.value, a.unit)
+Base.:-(a::RelativeQuantity{T, U}) where {T, U} = RelativeQuantity(-a.value, a.unit)
+
+# Scalar mul/div (Real to avoid ambiguity with unit-bearing types)
+Base.:*(a::Real, b::RelativeQuantity{T, U}) where {T, U} =
+    RelativeQuantity(a * b.value, b.unit)
+Base.:*(a::RelativeQuantity{T, U}, b::Real) where {T, U} =
+    RelativeQuantity(a.value * b, a.unit)
+Base.:/(a::RelativeQuantity{T, U}, b::Real) where {T, U} =
+    RelativeQuantity(a.value / b, a.unit)
+
+# Comparisons
+Base.:(==)(a::RelativeQuantity{T, U}, b::RelativeQuantity{S, U}) where {T, S, U} =
+    a.value == b.value
+Base.:(<)(a::RelativeQuantity{T, U}, b::RelativeQuantity{S, U}) where {T, S, U} =
+    a.value < b.value
+Base.:(<=)(a::RelativeQuantity{T, U}, b::RelativeQuantity{S, U}) where {T, S, U} =
+    a.value <= b.value
+Base.isless(a::RelativeQuantity{T, U}, b::RelativeQuantity{S, U}) where {T, S, U} =
+    isless(a.value, b.value)
+Base.isapprox(
+    a::RelativeQuantity{T, U},
+    b::RelativeQuantity{S, U};
+    kwargs...,
+) where {T, S, U} = isapprox(a.value, b.value; kwargs...)
+
+"""
+    ustrip(q::RelativeQuantity)
+
+Extract the numeric value from a `RelativeQuantity`.
+"""
+ustrip(q::RelativeQuantity) = q.value
+
+# Type conversions
+Base.convert(::Type{RelativeQuantity{T, U}}, q::RelativeQuantity{S, U}) where {T, S, U} =
+    RelativeQuantity(convert(T, q.value), q.unit)
+Base.promote_rule(
+    ::Type{RelativeQuantity{T, U}},
+    ::Type{RelativeQuantity{S, U}},
+) where {T, S, U} = RelativeQuantity{promote_type(T, S), U}
+
+# Display
+Base.show(io::IO, q::RelativeQuantity{T, DeviceBaseUnit}) where {T} =
+    print(io, q.value, " DU")
+Base.show(io::IO, q::RelativeQuantity{T, SystemBaseUnit}) where {T} =
+    print(io, q.value, " SU")
+Base.show(io::IO, ::DeviceBaseUnit) = print(io, "DU")
+Base.show(io::IO, ::SystemBaseUnit) = print(io, "SU")
+Base.show(io::IO, ::NaturalUnit) = print(io, "NU")
+
+Base.zero(::Type{RelativeQuantity{T, U}}) where {T, U} = RelativeQuantity(zero(T), U())
+Base.one(::Type{RelativeQuantity{T, U}}) where {T, U} = RelativeQuantity(one(T), U())
+
+"""
+    convert_cost_coefficient(value, U_from, U_to,
+                             system_base_power, device_base_power,
+                             exponent::Int = 1) → Float64
+
+Convert a cost coefficient (e.g. \$/MW for `exponent=1`, \$/MW² for
+`exponent=2`) between unit systems. The conversion ratio is the inverse of
+the corresponding power-value ratio raised to `exponent`, since if
+`obj = c · x_from` and `x_from = r · x_to`, then the equivalent coefficient
+under `x_to` is `c · r`.
+"""
+convert_cost_coefficient(
+    value::Float64,
+    U_from::AbstractUnitSystem,
+    U_to::AbstractUnitSystem,
+    system_base_power::Float64,
+    device_base_power::Float64,
+    exponent::Int = 1,
+) =
+    value * _cost_coeff_ratio(U_from, U_to, system_base_power, device_base_power)^exponent
+
+_cost_coeff_ratio(::SystemBaseUnit, ::SystemBaseUnit, _, _) = 1.0
+_cost_coeff_ratio(::DeviceBaseUnit, ::DeviceBaseUnit, _, _) = 1.0
+_cost_coeff_ratio(::NaturalUnit, ::NaturalUnit, _, _) = 1.0
+_cost_coeff_ratio(::DeviceBaseUnit, ::SystemBaseUnit, sb, db) = sb / db
+_cost_coeff_ratio(::SystemBaseUnit, ::DeviceBaseUnit, sb, db) = db / sb
+_cost_coeff_ratio(::NaturalUnit, ::SystemBaseUnit, sb, _) = sb
+_cost_coeff_ratio(::SystemBaseUnit, ::NaturalUnit, sb, _) = 1 / sb
+_cost_coeff_ratio(::NaturalUnit, ::DeviceBaseUnit, _, db) = db
+_cost_coeff_ratio(::DeviceBaseUnit, ::NaturalUnit, _, db) = 1 / db
+
+"""
+    display_units_arg(f, ::Type{T}) -> Union{AbstractRelativeUnit, Missing}
+
+Trait returning the units argument a getter `f` expects when called on a
+component of type `T` for display/tabular output, or `missing` if the getter
+takes no units argument. Keyed on both function and type because the same
+getter name can appear on both unit-bearing and non-unit-bearing structs
+(e.g. `get_b` on `Line` vs. `DynamicExponentialLoad`). Downstream packages
+set this per-struct (typically via the struct-generator template); consumers
+like `show_components` dispatch on the result to avoid runtime method
+introspection.
+"""
+display_units_arg(_, ::Type) = missing
diff --git a/src/time_series_interface.jl b/src/time_series_interface.jl
index 462f2b396..05ecfe040 100644
--- a/src/time_series_interface.jl
+++ b/src/time_series_interface.jl
@@ -950,7 +950,10 @@ function _make_time_array(owner, time_series, start_time, len, ignore_scaling_fa
         return ta
     end
 
-    return ta .* multiplier(owner)
+    # Scaling-factor multipliers (e.g. `get_max_active_power`) are unit-aware
+    # accessors from downstream packages; pass `SU` so the result is in the
+    # system base that consumers of the time series expect.
+    return ta .* multiplier(owner, SU)
 end
 
 """
diff --git a/src/utils/generate_structs.jl b/src/utils/generate_structs.jl
index 5dd0552d5..f29727fa6 100644
--- a/src/utils/generate_structs.jl
+++ b/src/utils/generate_structs.jl
@@ -63,13 +63,26 @@ end
 
 {{/has_null_values}}
 {{#accessors}}
+{{#needs_conversion}}
+{{#create_docstring}}\"\"\"Get [`{{struct_name}}`](@ref) `{{name}}`. The `units` argument is required (e.g. `SU`, `DU`, `MW`, or `Float64`).\"\"\"{{/create_docstring}}
+{{accessor}}(value::{{struct_name}}, units) = get_value(value, Val(:{{name}}), Val({{conversion_unit}}), units)
+InfrastructureSystems.display_units_arg(::typeof({{accessor}}), ::Type{ {{struct_name}} }) = InfrastructureSystems.SU
+{{/needs_conversion}}
+{{^needs_conversion}}
 {{#create_docstring}}\"\"\"Get [`{{struct_name}}`](@ref) `{{name}}`.\"\"\"{{/create_docstring}}
-{{accessor}}(value::{{struct_name}}) = {{#needs_conversion}}get_value(value, Val(:{{name}}), Val({{conversion_unit}})){{/needs_conversion}}{{^needs_conversion}}value.{{name}}{{/needs_conversion}}
+{{accessor}}(value::{{struct_name}}) = value.{{name}}
+{{/needs_conversion}}
 {{/accessors}}
 
 {{#setters}}
+{{#needs_conversion}}
 {{#create_docstring}}\"\"\"Set [`{{struct_name}}`](@ref) `{{name}}`.\"\"\"{{/create_docstring}}
-{{setter}}(value::{{struct_name}}, val) = value.{{name}} = {{#needs_conversion}}set_value(value, Val(:{{name}}), val, Val({{conversion_unit}})){{/needs_conversion}}{{^needs_conversion}}val{{/needs_conversion}}
+{{setter}}(value::{{struct_name}}, val) = value.{{name}} = set_value(value, Val(:{{name}}), val, Val({{conversion_unit}}))
+{{/needs_conversion}}
+{{^needs_conversion}}
+{{#create_docstring}}\"\"\"Set [`{{struct_name}}`](@ref) `{{name}}`.\"\"\"{{/create_docstring}}
+{{setter}}(value::{{struct_name}}, val) = value.{{name}} = val
+{{/needs_conversion}}
 {{/setters}}
 
 {{#custom_code}}
@@ -140,16 +153,32 @@ function generate_structs(directory, data::Vector; print_results = true)
             end
             accessor_name = accessor_module * "get_" * param["name"]
             setter_name = accessor_module * "set_" * param["name"] * "!"
-            push!(
-                accessors,
-                Dict(
-                    "name" => param["name"],
-                    "accessor" => accessor_name,
-                    "create_docstring" => create_docstring,
-                    "needs_conversion" => get(param, "needs_conversion", false),
-                    "conversion_unit" => get(param, "conversion_unit", "nothing"),
-                ),
-            )
+            conversion_unit = get(param, "conversion_unit", "nothing")
+            include_getter = !get(param, "exclude_getter", false)
+            if include_getter
+                push!(
+                    accessors,
+                    Dict(
+                        "name" => param["name"],
+                        "accessor" => accessor_name,
+                        "create_docstring" => create_docstring,
+                        "needs_conversion" => get(param, "needs_conversion", false),
+                        "conversion_unit" => conversion_unit,
+                    ),
+                )
+            else
+                internal_name = "_get_" * param["name"]
+                push!(
+                    accessors,
+                    Dict(
+                        "name" => param["name"],
+                        "accessor" => internal_name,
+                        "create_docstring" => false,
+                        "needs_conversion" => false,
+                        "conversion_unit" => "nothing",
+                    ),
+                )
+            end
             include_setter = !get(param, "exclude_setter", false)
             if include_setter
                 push!(
@@ -160,11 +189,14 @@ function generate_structs(directory, data::Vector; print_results = true)
                         "data_type" => param["data_type"],
                         "create_docstring" => create_docstring,
                         "needs_conversion" => get(param, "needs_conversion", false),
-                        "conversion_unit" => get(param, "conversion_unit", "nothing"),
+                        "conversion_unit" => conversion_unit,
                     ),
                 )
             end
             if field["name"] != "internal" && accessor_module == ""
+                # exclude_getter/exclude_setter mean "hand-written elsewhere" (e.g.
+                # unit-aware accessors with different signatures), not "nonexistent" —
+                # always export the public name.
                 push!(unique_accessor_functions, accessor_name)
                 push!(unique_setter_functions, setter_name)
             end
diff --git a/src/utils/print_pt_v2.jl b/src/utils/print_pt_v2.jl
index abf293a87..d8a3fd53e 100644
--- a/src/utils/print_pt_v2.jl
+++ b/src/utils/print_pt_v2.jl
@@ -148,7 +148,12 @@ function show_components(
                     val = summary(val)
                 elseif hasproperty(parent, getter_name)
                     getter_func = Base.getproperty(parent, getter_name)
-                    val = getter_func(component)
+                    arg = display_units_arg(getter_func, typeof(component))
+                    val = if ismissing(arg)
+                        getter_func(component)
+                    else
+                        getter_func(component, arg)
+                    end
                 end
                 data[i, j] = val
                 j += 1
diff --git a/src/utils/print_pt_v3.jl b/src/utils/print_pt_v3.jl
index 464c34bf4..e9d0ce8bf 100644
--- a/src/utils/print_pt_v3.jl
+++ b/src/utils/print_pt_v3.jl
@@ -129,7 +129,12 @@ function show_components(
                     val = summary(val)
                 elseif hasproperty(parent, getter_name)
                     getter_func = Base.getproperty(parent, getter_name)
-                    val = getter_func(component)
+                    arg = display_units_arg(getter_func, typeof(component))
+                    val = if ismissing(arg)
+                        getter_func(component)
+                    else
+                        getter_func(component, arg)
+                    end
                 end
                 data[i, j] = val
                 j += 1
diff --git a/src/utils/test.jl b/src/utils/test.jl
index 72e4bfba4..a14195daa 100644
--- a/src/utils/test.jl
+++ b/src/utils/test.jl
@@ -45,7 +45,11 @@ end
 get_internal(component::TestComponent) = component.internal
 get_internal(component::AdditionalTestComponent) = component.internal
 get_val(component::TestComponent) = component.val
+# 2-arg form so this getter can be used as a `scaling_factor_multiplier`
+# (which `_make_time_array` invokes with a units marker).
+get_val(component::TestComponent, _) = component.val
 get_val2(component::TestComponent) = component.val2
+get_val2(component::TestComponent, _) = component.val2
 supports_time_series(::TestComponent) = true
 supports_time_series(::AdditionalTestComponent) = true
 supports_time_series(::SimpleTestComponent) = false
diff --git a/src/value_curve.jl b/src/value_curve.jl
index 81efebbff..3b73246fc 100644
--- a/src/value_curve.jl
+++ b/src/value_curve.jl
@@ -272,12 +272,9 @@ end
 IncrementalCurve(data::AverageRateCurve) = IncrementalCurve(InputOutputCurve(data))
 
 # PRINTING
-"Whether there is a cost alias for the instance or type under consideration"
-is_cost_alias(::Union{ValueCurve, Type{<:ValueCurve}}) = false
-
-# For cost aliases, return the alias name; otherwise, return the type name without the parameter
-simple_type_name(curve::ValueCurve) =
-    string(is_cost_alias(curve) ? typeof(curve) : nameof(typeof(curve)))
+# typeof() can't recover const alias names, so we use nameof for non-aliases
+# and override in cost_aliases.jl for each alias.
+simple_type_name(curve::ValueCurve) = string(nameof(typeof(curve)))
 
 function Base.show(io::IO, ::MIME"text/plain", curve::InputOutputCurve)
     print(io, simple_type_name(curve))
diff --git a/test/test_cost_functions.jl b/test/test_cost_functions.jl
index 7a66b2612..22290eade 100644
--- a/test/test_cost_functions.jl
+++ b/test/test_cost_functions.jl
@@ -286,31 +286,35 @@ end
     @test zero(IS.FuelCurve) ==
           IS.FuelCurve(IS.InputOutputCurve(IS.LinearFunctionData(0.0, 0.0)), 0.0)
 
-    @test repr(cc) == sprint(show, cc) ==
-          "InfrastructureSystems.CostCurve{QuadraticCurve}(QuadraticCurve(1.0, 2.0, 3.0), InfrastructureSystems.UnitSystemModule.UnitSystem.NATURAL_UNITS = 2, LinearCurve(0.0, 0.0))"
-    @test repr(fc) == sprint(show, fc) ==
-          "InfrastructureSystems.FuelCurve{QuadraticCurve}(QuadraticCurve(1.0, 2.0, 3.0), InfrastructureSystems.UnitSystemModule.UnitSystem.NATURAL_UNITS = 2, 4.0, LinearCurve(0.0, 0.0), LinearCurve(0.0, 0.0))"
+    # repr and sprint(show, ...) must agree; the type parameter may or may not
+    # be module-qualified depending on what's in scope, so check key content.
+    @test repr(cc) == sprint(show, cc)
+    @test occursin("CostCurve", repr(cc))
+    @test occursin("QuadraticCurve(1.0, 2.0, 3.0)", repr(cc))
+    @test occursin("LinearCurve(0.0, 0.0)", repr(cc))
+    @test repr(fc) == sprint(show, fc)
+    @test occursin("FuelCurve", repr(fc))
+    @test occursin("QuadraticCurve(1.0, 2.0, 3.0)", repr(fc))
+    @test occursin("4.0", repr(fc))
     @test sprint(show, "text/plain", cc) ==
           sprint(show, "text/plain", cc; context = :compact => false) ==
-          "CostCurve:\n  value_curve: QuadraticCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x^2 + 2.0 x + 3.0\n  power_units: InfrastructureSystems.UnitSystemModule.UnitSystem.NATURAL_UNITS = 2\n  vom_cost: LinearCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 0.0 x + 0.0"
+          "CostCurve:\n  value_curve: QuadraticCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x^2 + 2.0 x + 3.0\n  vom_cost: LinearCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 0.0 x + 0.0\n  power_units: NU"
     @test sprint(show, "text/plain", fc) ==
           sprint(show, "text/plain", fc; context = :compact => false) ==
-          "FuelCurve:\n  value_curve: QuadraticCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x^2 + 2.0 x + 3.0\n  power_units: InfrastructureSystems.UnitSystemModule.UnitSystem.NATURAL_UNITS = 2\n  fuel_cost: 4.0\n  startup_fuel_offtake: LinearCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 0.0 x + 0.0\n  vom_cost: LinearCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 0.0 x + 0.0"
+          "FuelCurve:\n  value_curve: QuadraticCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x^2 + 2.0 x + 3.0\n  fuel_cost: 4.0\n  startup_fuel_offtake: LinearCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 0.0 x + 0.0\n  vom_cost: LinearCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 0.0 x + 0.0\n  power_units: NU"
     @test sprint(show, "text/plain", cc; context = :compact => true) ==
-          "CostCurve with power_units InfrastructureSystems.UnitSystemModule.UnitSystem.NATURAL_UNITS = 2, vom_cost LinearCurve(0.0, 0.0), and value_curve:\n  QuadraticCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x^2 + 2.0 x + 3.0"
+          "CostCurve with power_units NU, vom_cost LinearCurve(0.0, 0.0), and value_curve:\n  QuadraticCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x^2 + 2.0 x + 3.0"
     @test sprint(show, "text/plain", fc; context = :compact => true) ==
-          "FuelCurve with power_units InfrastructureSystems.UnitSystemModule.UnitSystem.NATURAL_UNITS = 2, fuel_cost 4.0, startup_fuel_offtake LinearCurve(0.0, 0.0), vom_cost LinearCurve(0.0, 0.0), and value_curve:\n  QuadraticCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x^2 + 2.0 x + 3.0"
+          "FuelCurve with power_units NU, fuel_cost 4.0, startup_fuel_offtake LinearCurve(0.0, 0.0), vom_cost LinearCurve(0.0, 0.0), and value_curve:\n  QuadraticCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x^2 + 2.0 x + 3.0"
 
-    @test IS.get_power_units(cc) == IS.UnitSystem.NATURAL_UNITS
-    @test IS.get_power_units(fc) == IS.UnitSystem.NATURAL_UNITS
+    @test IS.get_power_units(cc) == IS.NaturalUnit()
+    @test IS.get_power_units(fc) == IS.NaturalUnit()
     @test IS.get_power_units(
-        IS.CostCurve(zero(IS.InputOutputCurve), IS.UnitSystem.SYSTEM_BASE),
-    ) ==
-          IS.UnitSystem.SYSTEM_BASE
+        IS.CostCurve(zero(IS.InputOutputCurve), IS.SystemBaseUnit()),
+    ) == IS.SystemBaseUnit()
     @test IS.get_power_units(
-        IS.FuelCurve(zero(IS.InputOutputCurve), IS.UnitSystem.DEVICE_BASE, 1.0),
-    ) ==
-          IS.UnitSystem.DEVICE_BASE
+        IS.FuelCurve(zero(IS.InputOutputCurve), IS.DeviceBaseUnit(), 1.0),
+    ) == IS.DeviceBaseUnit()
 
     @test IS.get_vom_cost(cc) == IS.LinearCurve(0.0)
     @test IS.get_vom_cost(fc) == IS.LinearCurve(0.0)
diff --git a/test/test_make_convex.jl b/test/test_make_convex.jl
index 68016cd09..e7361293e 100644
--- a/test/test_make_convex.jl
+++ b/test/test_make_convex.jl
@@ -198,14 +198,14 @@
         vom_cost = IS.LinearCurve(0.5)
         cost_curve_convex = IS.CostCurve(;
             value_curve = ioc_convex,
-            power_units = IS.UnitSystem.NATURAL_UNITS,
+            power_units = IS.NaturalUnit(),
             vom_cost = vom_cost,
         )
 
         result = IS.increasing_curve_convex_approximation(cost_curve_convex)
         @test result isa IS.CostCurve
         @test IS.is_convex(result)
-        @test IS.get_power_units(result) == IS.UnitSystem.NATURAL_UNITS
+        @test IS.get_power_units(result) == IS.NaturalUnit()
         @test IS.get_vom_cost(result) == vom_cost
 
         # Non-convex CostCurve - should be convexified
@@ -218,14 +218,14 @@
         vom_cost_2 = IS.LinearCurve(1.0)
         cost_curve_concave = IS.CostCurve(;
             value_curve = ioc_concave,
-            power_units = IS.UnitSystem.SYSTEM_BASE,
+            power_units = IS.SystemBaseUnit(),
             vom_cost = vom_cost_2,
         )
 
         result_concave = IS.increasing_curve_convex_approximation(cost_curve_concave)
         @test result_concave isa IS.CostCurve
         @test IS.is_convex(result_concave)
-        @test IS.get_power_units(result_concave) == IS.UnitSystem.SYSTEM_BASE
+        @test IS.get_power_units(result_concave) == IS.SystemBaseUnit()
         @test IS.get_vom_cost(result_concave) == vom_cost_2
 
         # Invalid CostCurve (not strictly increasing) - should throw error
@@ -256,7 +256,7 @@
         vom_cost = IS.LinearCurve(0.5)
         fuel_curve_convex = IS.FuelCurve(;
             value_curve = ioc_convex,
-            power_units = IS.UnitSystem.NATURAL_UNITS,
+            power_units = IS.NaturalUnit(),
             fuel_cost = 25.0,
             vom_cost = vom_cost,
         )
@@ -264,7 +264,7 @@
         result = IS.increasing_curve_convex_approximation(fuel_curve_convex)
         @test result isa IS.FuelCurve
         @test IS.is_convex(result)
-        @test IS.get_power_units(result) == IS.UnitSystem.NATURAL_UNITS
+        @test IS.get_power_units(result) == IS.NaturalUnit()
         @test result.fuel_cost == 25.0
         @test IS.get_vom_cost(result) == vom_cost
 
@@ -278,7 +278,7 @@
         vom_cost_2 = IS.LinearCurve(1.0)
         fuel_curve_concave = IS.FuelCurve(;
             value_curve = ioc_concave,
-            power_units = IS.UnitSystem.SYSTEM_BASE,
+            power_units = IS.SystemBaseUnit(),
             fuel_cost = 30.0,
             vom_cost = vom_cost_2,
         )
@@ -286,7 +286,7 @@
         result_concave = IS.increasing_curve_convex_approximation(fuel_curve_concave)
         @test result_concave isa IS.FuelCurve
         @test IS.is_convex(result_concave)
-        @test IS.get_power_units(result_concave) == IS.UnitSystem.SYSTEM_BASE
+        @test IS.get_power_units(result_concave) == IS.SystemBaseUnit()
         @test result_concave.fuel_cost == 30.0
         @test IS.get_vom_cost(result_concave) == vom_cost_2
 
@@ -295,7 +295,7 @@
         inc_convex = IS.IncrementalCurve(psd_convex, 0.0)
         fuel_curve_inc = IS.FuelCurve(;
             value_curve = inc_convex,
-            power_units = IS.UnitSystem.NATURAL_UNITS,
+            power_units = IS.NaturalUnit(),
             fuel_cost = 20.0,
         )
 
@@ -315,7 +315,7 @@
             ioc_invalid = IS.InputOutputCurve(pld_invalid)
             fuel_curve_invalid = IS.FuelCurve(;
                 value_curve = ioc_invalid,
-                power_units = IS.UnitSystem.NATURAL_UNITS,
+                power_units = IS.NaturalUnit(),
                 fuel_cost = 25.0,
             )
 
diff --git a/test/test_relative_units.jl b/test/test_relative_units.jl
new file mode 100644
index 000000000..4cae6bfff
--- /dev/null
+++ b/test/test_relative_units.jl
@@ -0,0 +1,79 @@
+@testset "RelativeQuantity construction and arithmetic" begin
+    a = 0.6 * IS.DU
+    b = 0.4 * IS.DU
+    @test a isa IS.RelativeQuantity{Float64, IS.DeviceBaseUnit}
+    @test IS.ustrip(a + b) ≈ 1.0
+    @test IS.ustrip(a - b) ≈ 0.2
+    @test IS.ustrip(-a) ≈ -0.6
+    # scalar multiplication dispatches differently on each side
+    @test IS.ustrip(2.0 * a) ≈ 1.2
+    @test IS.ustrip(a * 2.0) ≈ 1.2
+    @test IS.ustrip(a / 2.0) ≈ 0.3
+end
+
+@testset "RelativeQuantity comparisons" begin
+    @test 0.6 * IS.DU < 0.7 * IS.DU
+    @test 0.6 * IS.DU <= 0.6 * IS.DU
+    @test isapprox(0.6 * IS.DU, 0.60000001 * IS.DU; atol = 1e-6)
+    @test isless(0.6 * IS.DU, 0.7 * IS.DU)
+end
+
+@testset "DU and SU cannot be mixed" begin
+    @test_throws Exception 0.6 * IS.DU + 0.4 * IS.SU
+    @test_throws Exception 0.6 * IS.DU == 0.4 * IS.SU
+end
+
+@testset "RelativeQuantity zero and one" begin
+    @test zero(IS.RelativeQuantity{Float64, IS.DeviceBaseUnit}) == 0.0 * IS.DU
+    @test one(IS.RelativeQuantity{Float64, IS.DeviceBaseUnit}) == 1.0 * IS.DU
+end
+
+@testset "RelativeQuantity display" begin
+    @test sprint(show, 0.6 * IS.DU) == "0.6 DU"
+    @test sprint(show, 0.3 * IS.SU) == "0.3 SU"
+    @test sprint(show, IS.DU) == "DU"
+    @test sprint(show, IS.SU) == "SU"
+    @test sprint(show, IS.NU) == "NU"
+end
+
+@testset "convert_cost_coefficient" begin
+    sb, db = 100.0, 50.0
+    @testset "identity (same unit system)" begin
+        for U in (IS.SU, IS.DU, IS.NU)
+            @test IS.convert_cost_coefficient(2.5, U, U, sb, db) == 2.5
+            @test IS.convert_cost_coefficient(2.5, U, U, sb, db, 2) == 2.5
+        end
+    end
+
+    @testset "DU ↔ SU (linear)" begin
+        @test IS.convert_cost_coefficient(2.0, IS.DU, IS.SU, sb, db) ≈ 2.0 * sb / db
+        @test IS.convert_cost_coefficient(2.0, IS.SU, IS.DU, sb, db) ≈ 2.0 * db / sb
+    end
+
+    @testset "NU ↔ {SU, DU} (linear)" begin
+        @test IS.convert_cost_coefficient(2.0, IS.NU, IS.SU, sb, db) ≈ 2.0 * sb
+        @test IS.convert_cost_coefficient(2.0, IS.SU, IS.NU, sb, db) ≈ 2.0 / sb
+        @test IS.convert_cost_coefficient(2.0, IS.NU, IS.DU, sb, db) ≈ 2.0 * db
+        @test IS.convert_cost_coefficient(2.0, IS.DU, IS.NU, sb, db) ≈ 2.0 / db
+    end
+
+    @testset "exponent (quadratic)" begin
+        @test IS.convert_cost_coefficient(2.0, IS.DU, IS.SU, sb, db, 2) ≈ 2.0 * (sb / db)^2
+        @test IS.convert_cost_coefficient(2.0, IS.NU, IS.SU, sb, db, 2) ≈ 2.0 * sb^2
+    end
+
+    @testset "round-trip is identity (linear and quadratic)" begin
+        for (Ua, Ub) in ((IS.DU, IS.SU), (IS.NU, IS.SU), (IS.NU, IS.DU))
+            for k in (1, 2)
+                forward = IS.convert_cost_coefficient(2.0, Ua, Ub, sb, db, k)
+                back = IS.convert_cost_coefficient(forward, Ub, Ua, sb, db, k)
+                @test back ≈ 2.0
+            end
+        end
+    end
+
+    @testset "negative exponent inverts linear ratio (used for piecewise x-coords)" begin
+        @test IS.convert_cost_coefficient(2.0, IS.DU, IS.SU, sb, db, -1) ≈ 2.0 * db / sb
+        @test IS.convert_cost_coefficient(2.0, IS.NU, IS.SU, sb, db, -1) ≈ 2.0 / sb
+    end
+end
diff --git a/test/test_time_series_function_data.jl b/test/test_time_series_function_data.jl
index cc9987eef..cb724e018 100644
--- a/test/test_time_series_function_data.jl
+++ b/test/test_time_series_function_data.jl
@@ -386,9 +386,9 @@ end
             IS.TimeSeriesLinearFunctionData(forecast_key),
         )
         # CostCurve preserves value_curve and accepts power_units
-        cc = IS.CostCurve(ts_io, IS.UnitSystem.SYSTEM_BASE)
+        cc = IS.CostCurve(ts_io, IS.SystemBaseUnit())
         @test IS.get_value_curve(cc) === ts_io
-        @test IS.get_power_units(cc) == IS.UnitSystem.SYSTEM_BASE
+        @test IS.get_power_units(cc) == IS.SystemBaseUnit()
 
         # FuelCurve preserves fuel_cost
         fc = IS.FuelCurve(ts_io, 5.0)