diff --git a/Project.toml b/Project.toml
index 0a828b7..0c2127c 100644
--- a/Project.toml
+++ b/Project.toml
@@ -26,7 +26,7 @@ PowerFlows = "94fada2c-0ca5-4b90-a1fb-4bc5b59ccfc7"
 [sources]
 InfrastructureSystems = {rev = "IS4", url = "https://github.com/Sienna-Platform/InfrastructureSystems.jl"}
 PowerSystems = {rev = "psy6", url = "https://github.com/Sienna-Platform/PowerSystems.jl"}
-InfrastructureOptimizationModels = {rev = "main", url = "https://github.com/Sienna-Platform/InfrastructureOptimizationModels.jl"}
+InfrastructureOptimizationModels = {rev = "ac/hvdc-vsc", url = "https://github.com/Sienna-Platform/InfrastructureOptimizationModels.jl"}
 
 [extensions]
 PowerFlowsExt = "PowerFlows"
@@ -34,7 +34,6 @@ PowerFlowsExt = "PowerFlows"
 [compat]
 Dates = "1"
 DocStringExtensions = "~0.8, ~0.9"
-InfrastructureOptimizationModels = "0.1.0"
 InfrastructureSystems = "3"
 InteractiveUtils = "1.11.0"
 JuMP = "^1.28"
diff --git a/src/PowerOperationsModels.jl b/src/PowerOperationsModels.jl
index 738cc21..6da666b 100644
--- a/src/PowerOperationsModels.jl
+++ b/src/PowerOperationsModels.jl
@@ -213,6 +213,8 @@ include("common_models/add_to_expression.jl")
 include("common_models/add_parameters.jl")
 include("common_models/make_system_expressions.jl")
 include("common_models/reserve_range_constraints.jl")
+include("common_models/quadratic_converter_loss.jl")
+include("common_models/network_conditional.jl")
 
 # Market bid cost plumbing (PSY orchestration moved out of IOM). Must be included
 # before device-specific files that reference MBC_TYPES / IEC_TYPES.
@@ -506,11 +508,15 @@ export PostContingencyActivePowerReserveDeploymentVariable
 export DCVoltage
 export DCLineCurrent
 export ConverterCurrent
-export ConverterPositiveCurrent
-export ConverterNegativeCurrent
-export ConverterCurrentDirection
+export PositiveCurrent
+export NegativeCurrent
+export CurrentDirection
 export HVDCFlowDirectionVariable
 export HVDCLosses
+export HVDCFromDCVoltage
+export HVDCToDCVoltage
+export HVDCReactivePowerFromVariable
+export HVDCReactivePowerToVariable
 
 # Load Variables
 export ShiftUpActivePowerVariable
@@ -758,12 +764,14 @@ export HVDCTwoTerminalLossless
 export HVDCTwoTerminalDispatch
 export HVDCTwoTerminalPiecewiseLoss
 export HVDCTwoTerminalLCC
+export HVDCTwoTerminalVSCNLP
+export HVDCTwoTerminalVSCLP
 
 # Converter Formulations
 export LosslessConverter
 export LinearLossConverter
 export AbstractQuadraticLossConverter
-export Bin2QuadraticLossConverter
+export MILPQuadraticLossConverter
 export QuadraticLossConverter
 
 # DC Line Formulations
diff --git a/src/ac_transmission_models/branch_constructor.jl b/src/ac_transmission_models/branch_constructor.jl
index 568b072..fcc96db 100644
--- a/src/ac_transmission_models/branch_constructor.jl
+++ b/src/ac_transmission_models/branch_constructor.jl
@@ -1670,3 +1670,161 @@ function construct_device!(
     add_feedforward_constraints!(container, device_model, devices)
     return
 end
+
+############################################################################
+####################### Two-Terminal VSC HVDC Construct ####################
+############################################################################
+
+# Quadratic / bilinear approximation traits — same scheme used by the MT
+# converter formulations.
+_quad_config(::Type{HVDCTwoTerminalVSCNLP}) = IOM.NoQuadApproxConfig()
+_quad_config(::Type{HVDCTwoTerminalVSCLP}) =
+    IOM.SolverSOS2QuadConfig(DEFAULT_INTERPOLATION_LENGTH)
+_bilinear_config(::Type{HVDCTwoTerminalVSCNLP}) = IOM.NoBilinearApproxConfig()
+_bilinear_config(::Type{HVDCTwoTerminalVSCLP}) =
+    IOM.Bin2Config(IOM.SolverSOS2QuadConfig(DEFAULT_INTERPOLATION_LENGTH))
+
+function construct_device!(
+    container::OptimizationContainer,
+    sys::PSY.System,
+    ::ArgumentConstructStage,
+    device_model::DeviceModel{PSY.TwoTerminalVSCLine, F},
+    network_model::NetworkModel{<:AbstractPowerModel},
+) where {F <: AbstractTwoTerminalVSCFormulation}
+    devices = get_available_components(device_model, sys)
+
+    add_variables!(container, FlowActivePowerFromToVariable, devices, F)
+    add_variables!(container, FlowActivePowerToFromVariable, devices, F)
+    add_variables!(container, DCLineCurrentFlowVariable, devices, F)
+    add_variables!(container, HVDCFromDCVoltage, devices, F)
+    add_variables!(container, HVDCToDCVoltage, devices, F)
+
+    _maybe_add_reactive_power_variables!(
+        container, devices, device_model, network_model,
+        (HVDCReactivePowerFromVariable, HVDCReactivePowerToVariable),
+    )
+
+    # use_linear_loss = true adds a binary direction variable (CurrentDirection).
+    # Both HVDCTwoTerminalVSCNLP and HVDCTwoTerminalVSCLP dispatch through here:
+    # on the NLP formulation this pushes the model from a smooth NLP to MINLP
+    # (caller must supply a MINLP-capable solver); on the LP formulation the
+    # SOS2 PWL approximations already make the model MILP, so the extra binary
+    # is free.
+    if get_attribute(device_model, "use_linear_loss")
+        _add_abs_value_decomposition_variables!(
+            container, devices, device_model, network_model,
+        )
+    end
+
+    add_to_expression!(
+        container, ActivePowerBalance, FlowActivePowerFromToVariable,
+        devices, device_model, network_model,
+    )
+    add_to_expression!(
+        container, ActivePowerBalance, FlowActivePowerToFromVariable,
+        devices, device_model, network_model,
+    )
+
+    add_feedforward_arguments!(container, device_model, devices)
+    return
+end
+
+function construct_device!(
+    container::OptimizationContainer,
+    sys::PSY.System,
+    ::ModelConstructStage,
+    device_model::DeviceModel{PSY.TwoTerminalVSCLine, F},
+    network_model::NetworkModel{<:AbstractPowerModel},
+) where {F <: AbstractTwoTerminalVSCFormulation}
+    devices = get_available_components(device_model, sys)
+    time_steps = get_time_steps(container)
+    line_names = [PSY.get_name(d) for d in devices]
+
+    v_f_var = get_variable(container, HVDCFromDCVoltage, PSY.TwoTerminalVSCLine)
+    v_t_var = get_variable(container, HVDCToDCVoltage, PSY.TwoTerminalVSCLine)
+    i_var = get_variable(container, DCLineCurrentFlowVariable, PSY.TwoTerminalVSCLine)
+
+    v_f_bounds = PSY.get_voltage_limits_from.(devices)
+    v_t_bounds = PSY.get_voltage_limits_to.(devices)
+    i_bounds = [
+        (min = -_linear_loss_i_max(d), max = _linear_loss_i_max(d)) for d in devices
+    ]
+
+    quad_cfg, bilin_cfg = _quad_config(F), _bilinear_config(F)
+
+    v_f_sq_expr = IOM._add_quadratic_approx!(
+        quad_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps, v_f_var, v_f_bounds, "v_f_sq",
+    )
+    v_t_sq_expr = IOM._add_quadratic_approx!(
+        quad_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps, v_t_var, v_t_bounds, "v_t_sq",
+    )
+    i_sq_expr = IOM._add_quadratic_approx!(
+        quad_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps, i_var, i_bounds, "i_sq",
+    )
+
+    IOM._add_bilinear_approx!(
+        bilin_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps,
+        v_f_sq_expr, i_sq_expr, v_f_var, i_var,
+        v_f_bounds, i_bounds, "vi_ft",
+    )
+    IOM._add_bilinear_approx!(
+        bilin_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps,
+        v_t_sq_expr, i_sq_expr, v_t_var, i_var,
+        v_t_bounds, i_bounds, "vi_tf",
+    )
+
+    _register_pq_sq_expressions!(
+        container, devices, line_names, time_steps, device_model,
+        network_model,
+    )
+
+    if get_attribute(device_model, "use_linear_loss")
+        _add_abs_value_decomposition_constraints!(
+            container, devices, device_model, network_model,
+            DCLineCurrentFlowVariable,
+        )
+    end
+
+    add_constraints!(
+        container, HVDCCableOhmsLawConstraint, devices, device_model, network_model,
+    )
+    add_constraints!(
+        container, HVDCVSCConverterPowerConstraint, devices, device_model, network_model,
+    )
+    _maybe_add_reactive_power_constraints!(
+        container, devices, device_model, network_model,
+        HVDCVSCApparentPowerLimitConstraint,
+    )
+
+    add_constraint_dual!(container, sys, device_model)
+    add_feedforward_constraints!(container, device_model, devices)
+    return
+end
+
+# AreaBalancePowerModel warning (consistent with other two-terminal formulations).
+function construct_device!(
+    ::OptimizationContainer,
+    ::PSY.System,
+    ::ArgumentConstructStage,
+    ::DeviceModel{PSY.TwoTerminalVSCLine, <:AbstractTwoTerminalVSCFormulation},
+    ::NetworkModel{AreaBalancePowerModel},
+)
+    @warn "AreaBalancePowerModel doesn't model individual line flows for PSY.TwoTerminalVSCLine. Arguments not built"
+    return
+end
+
+function construct_device!(
+    ::OptimizationContainer,
+    ::PSY.System,
+    ::ModelConstructStage,
+    ::DeviceModel{PSY.TwoTerminalVSCLine, <:AbstractTwoTerminalVSCFormulation},
+    ::NetworkModel{AreaBalancePowerModel},
+)
+    @warn "AreaBalancePowerModel doesn't model individual line flows for PSY.TwoTerminalVSCLine. Model not built"
+    return
+end
diff --git a/src/common_models/add_to_expression.jl b/src/common_models/add_to_expression.jl
index e8e9c7e..576f056 100644
--- a/src/common_models/add_to_expression.jl
+++ b/src/common_models/add_to_expression.jl
@@ -718,6 +718,45 @@ function add_to_expression!(
     return
 end
 
+# A VSC terminal can inject or consume Q freely, so the variable enters
+# `ReactivePowerBalance` as a signed injection (+1.0) rather than a load (−1.0).
+# Side selection picks the from- or to-terminal bus via dispatch on the
+# variable type, so the body is written once.
+_vsc_q_terminal_bus(d, ::Type{HVDCReactivePowerFromVariable}) = PSY.get_arc(d).from
+_vsc_q_terminal_bus(d, ::Type{HVDCReactivePowerToVariable}) = PSY.get_arc(d).to
+
+function add_to_expression!(
+    container::OptimizationContainer,
+    ::Type{T},
+    ::Type{U},
+    devices::IS.FlattenIteratorWrapper{V},
+    ::DeviceModel{V, W},
+    network_model::NetworkModel{X},
+) where {
+    T <: ReactivePowerBalance,
+    U <: Union{HVDCReactivePowerFromVariable, HVDCReactivePowerToVariable},
+    V <: PSY.TwoTerminalVSCLine,
+    W <: AbstractTwoTerminalVSCFormulation,
+    X <: ACPPowerModel,
+}
+    var = get_variable(container, U, V)
+    nodal_expr = get_expression(container, T, PSY.ACBus)
+    network_reduction = get_network_reduction(network_model)
+    time_steps = get_time_steps(container)
+    for d in devices
+        name = PSY.get_name(d)
+        bus_no = PNM.get_mapped_bus_number(
+            network_reduction, _vsc_q_terminal_bus(d, U),
+        )
+        for t in time_steps
+            add_proportional_to_jump_expression!(
+                nodal_expr[bus_no, t], var[name, t], 1.0,
+            )
+        end
+    end
+    return
+end
+
 """
 PWL implementation to add FromTo branch variables to SystemBalanceExpressions
 """
diff --git a/src/common_models/network_conditional.jl b/src/common_models/network_conditional.jl
new file mode 100644
index 0000000..e11a5c6
--- /dev/null
+++ b/src/common_models/network_conditional.jl
@@ -0,0 +1,59 @@
+# Helpers for variables/constraints that should only appear when the network
+# model actually represents the relevant physical quantity (e.g. reactive
+# power on AC networks). Each helper has a no-op method that dispatches on
+# `NetworkModel{<:AbstractActivePowerModel}`; Julia's method resolution picks
+# the more specific no-op over the AC body for DC formulations.
+
+"""
+Add reactive-power variables for a device and register them in the system's
+`ReactivePowerBalance` expression. `var_types` is a tuple/iterable of
+`VariableType` subtypes; each is added via `add_variables!` and then linked
+into `ReactivePowerBalance` via `add_to_expression!`. The caller's
+device-specific `add_to_expression!` methods are responsible for the actual
+bus mapping and sign convention.
+"""
+function _maybe_add_reactive_power_variables!(
+    container::OptimizationContainer,
+    devices,
+    model::DeviceModel{D, F},
+    network_model::NetworkModel{<:AbstractPowerModel},
+    var_types,
+) where {D <: PSY.Device, F}
+    for V in var_types
+        add_variables!(container, V, devices, F)
+        add_to_expression!(
+            container, ReactivePowerBalance, V, devices, model, network_model,
+        )
+    end
+    return
+end
+
+_maybe_add_reactive_power_variables!(
+    ::OptimizationContainer,
+    _devices,
+    ::DeviceModel{D, F},
+    ::NetworkModel{<:AbstractActivePowerModel},
+    _var_types,
+) where {D <: PSY.Device, F} = nothing
+
+"""
+Add a reactive-power-related constraint for a device on AC networks.
+"""
+function _maybe_add_reactive_power_constraints!(
+    container::OptimizationContainer,
+    devices,
+    model::DeviceModel{D, F},
+    network_model::NetworkModel{<:AbstractPowerModel},
+    constraint_type::Type{<:ConstraintType},
+) where {D <: PSY.Device, F}
+    add_constraints!(container, constraint_type, devices, model, network_model)
+    return
+end
+
+_maybe_add_reactive_power_constraints!(
+    ::OptimizationContainer,
+    _devices,
+    ::DeviceModel{D, F},
+    ::NetworkModel{<:AbstractActivePowerModel},
+    ::Type{<:ConstraintType},
+) where {D <: PSY.Device, F} = nothing
diff --git a/src/common_models/quadratic_converter_loss.jl b/src/common_models/quadratic_converter_loss.jl
new file mode 100644
index 0000000..9ad09d7
--- /dev/null
+++ b/src/common_models/quadratic_converter_loss.jl
@@ -0,0 +1,131 @@
+# Shared helpers for quadratic / two-term converter losses
+#   loss(I) = a * I^2 + b * |I| + c
+# Used by multi-terminal InterconnectingConverter formulations
+# (MILPQuadraticLossConverter, QuadraticLossConverter) and two-terminal
+# HVDCTwoTerminalVSCNLP formulations.
+
+#########################################
+######## Loss-curve introspection #######
+#########################################
+
+_get_quadratic_term(loss_fn::PSY.QuadraticCurve) = PSY.get_quadratic_term(loss_fn)
+_get_quadratic_term(loss_fn) = 0.0
+
+_has_linear_loss(d::PSY.InterconnectingConverter) =
+    !iszero(PSY.get_proportional_term(PSY.get_loss_function(d)))
+_has_linear_loss(d::PSY.TwoTerminalVSCLine) =
+    !iszero(PSY.get_proportional_term(PSY.get_converter_loss_from(d))) ||
+    !iszero(PSY.get_proportional_term(PSY.get_converter_loss_to(d)))
+
+# `filter` has no method for FlattenIteratorWrapper, so use a comprehension
+# that works on any iterable.
+_devices_with_linear_loss(devices) = [d for d in devices if _has_linear_loss(d)]
+
+#########################################
+######## Loss expression builder ########
+#########################################
+
+# `_loss_seed` picks the right JuMP expression flavor for the I^2 term up
+# front: QuadExpr when `i_sq_t` is the exact bilinear i*i (NLP path,
+# NoQuadApproxConfig), AffExpr when it's the SOS2 PWL surrogate (MIP path).
+# With `a == 0` the I^2 term drops, so we degrade to AffExpr in either case.
+_loss_seed(c::Float64, ::JuMP.QuadExpr) = JuMP.QuadExpr(JuMP.AffExpr(c))
+_loss_seed(c::Float64, ::JuMP.AffExpr) = JuMP.AffExpr(c)
+
+function _quadratic_converter_loss_expr(
+    a::Float64, b::Float64, c::Float64,
+    i_sq_t, i_pos_t, i_neg_t;
+    use_linear_loss::Bool,
+)
+    expr = iszero(a) ? JuMP.AffExpr(c) : _loss_seed(c, i_sq_t)
+    iszero(a) || JuMP.add_to_expression!(expr, a, i_sq_t)
+    if use_linear_loss && !iszero(b)
+        JuMP.add_to_expression!(expr, b, i_pos_t)
+        JuMP.add_to_expression!(expr, b, i_neg_t)
+    end
+    return expr
+end
+
+#########################################
+####### Abs-value decomposition #########
+#########################################
+
+# Split into two stages so the variables are added in ArgumentConstructStage
+# and the constraints in ModelConstructStage, matching the IOM construction
+# convention. Both no-op when no device has a nonzero proportional loss term;
+# the variables-stage helper additionally emits a configuration warning.
+function _add_abs_value_decomposition_variables!(
+    container::OptimizationContainer,
+    devices,
+    ::DeviceModel{D, F},
+    ::NetworkModel{<:AbstractPowerModel},
+) where {D <: PSY.Device, F}
+    ll_devices = _devices_with_linear_loss(devices)
+    if isempty(ll_devices)
+        @warn "use_linear_loss is enabled but no $(D) has a nonzero proportional loss term; no linear-loss variables/constraints will be added."
+        return
+    end
+    add_variables!(container, PositiveCurrent, ll_devices, F)
+    add_variables!(container, NegativeCurrent, ll_devices, F)
+    add_variables!(container, CurrentDirection, ll_devices, F)
+    return
+end
+
+# I_max comes from different PSY accessors depending on the device, so we
+# dispatch on the device type rather than asking the caller to thread a getter
+# through. Two terminals on a VSC line share the same DC current variable, so
+# we take the binding (min) rating; MT converters expose a single getter.
+_linear_loss_i_max(d::PSY.TwoTerminalVSCLine) =
+    min(PSY.get_max_dc_current_from(d), PSY.get_max_dc_current_to(d))
+_linear_loss_i_max(d::PSY.InterconnectingConverter) = PSY.get_max_dc_current(d)
+
+function _add_abs_value_decomposition_constraints!(
+    container::OptimizationContainer,
+    devices,
+    ::DeviceModel{D, F},
+    ::NetworkModel{<:AbstractPowerModel},
+    parent_var_type::Type{<:VariableType},
+) where {D <: PSY.Device, F}
+    ll_devices = _devices_with_linear_loss(devices)
+    isempty(ll_devices) && return
+
+    time_steps = get_time_steps(container)
+    names = [PSY.get_name(d) for d in ll_devices]
+    jump_model = get_jump_model(container)
+    i_var = get_variable(container, parent_var_type, D)
+    i_pos_var = get_variable(container, PositiveCurrent, D)
+    i_neg_var = get_variable(container, NegativeCurrent, D)
+    i_dir_var = get_variable(container, CurrentDirection, D)
+
+    abs_val_const = add_constraints_container!(
+        container, CurrentAbsoluteValueConstraint, D, names, time_steps,
+    )
+    pos_ub_const = add_constraints_container!(
+        container, CurrentAbsoluteValueConstraint, D, names, time_steps;
+        meta = "pos_ub",
+    )
+    neg_ub_const = add_constraints_container!(
+        container, CurrentAbsoluteValueConstraint, D, names, time_steps;
+        meta = "neg_ub",
+    )
+
+    for d in ll_devices
+        name = PSY.get_name(d)
+        i_max = _linear_loss_i_max(d)
+        for t in time_steps
+            abs_val_const[name, t] = JuMP.@constraint(
+                jump_model,
+                i_var[name, t] == i_pos_var[name, t] - i_neg_var[name, t],
+            )
+            pos_ub_const[name, t] = JuMP.@constraint(
+                jump_model,
+                i_pos_var[name, t] <= i_max * i_dir_var[name, t],
+            )
+            neg_ub_const[name, t] = JuMP.@constraint(
+                jump_model,
+                i_neg_var[name, t] <= i_max * (1 - i_dir_var[name, t]),
+            )
+        end
+    end
+    return
+end
diff --git a/src/core/constraints.jl b/src/core/constraints.jl
index 634e5c2..1eb459d 100644
--- a/src/core/constraints.jl
+++ b/src/core/constraints.jl
@@ -419,6 +419,42 @@ v_d^i = v_d^r - R_d I_d
 """
 struct HVDCTransmissionDCLineConstraint <: ConstraintType end
 
+"""
+Per-terminal converter power-balance constraint for two-terminal VSC HVDC:
+
+```math
+\\begin{aligned}
+p_{ft} &= v_f \\cdot I + (a_f I^2 + b_f |I| + c_f) \\\\
+p_{tf} &= -v_t \\cdot I + (a_t I^2 + b_t |I| + c_t)
+\\end{aligned}
+```
+"""
+struct HVDCVSCConverterPowerConstraint <: ConstraintType end
+
+"""
+Cable Ohm's law for a two-terminal HVDC link with explicit DC resistance:
+
+```math
+v_f - v_t = (1/g) \\cdot I
+```
+"""
+struct HVDCCableOhmsLawConstraint <: ConstraintType end
+
+"""
+Apparent-power limit at each terminal of a two-terminal VSC HVDC, added only on
+AC networks. Enforces ``|S_k| \\le S_k^{\\max}`` for ``k \\in \\{f, t\\}`` via one
+of two formulation-specific shapes:
+
+- `HVDCTwoTerminalVSCNLP` (NLP): exact disk ``p_k^2 + q_k^2 \\le (S_k^{\\max})^2``.
+- `HVDCTwoTerminalVSCLP` (LP): linear outer-approximation. The axis-aligned box
+  ``|p_k|, |q_k| \\le S_k^{\\max}`` is added unconditionally; the four diagonal
+  half-planes ``|p_k| \\pm q_k \\le S_k^{\\max}\\sqrt{2}`` are added when the
+  device-model attribute `use_octagon` (default `true`) is on, in which case
+  the intersection is a regular octagon circumscribing the disk
+  (loose by at most ≈8.2% in area).
+"""
+struct HVDCVSCApparentPowerLimitConstraint <: ConstraintType end
+
 abstract type PowerVariableLimitsConstraint <: ConstraintType end
 """
 Struct to create the constraint to limit active power input expressions.
diff --git a/src/core/formulations.jl b/src/core/formulations.jl
index c37aec3..32201ae 100644
--- a/src/core/formulations.jl
+++ b/src/core/formulations.jl
@@ -184,8 +184,33 @@ Branch type to represent non-linear LCC (line commutated converter) model on two
 """
 struct HVDCTwoTerminalLCC <: AbstractTwoTerminalDCLineFormulation end
 
-# Not Implemented
-# struct VoltageSourceDC <: AbstractTwoTerminalDCLineFormulation end
+"""
+Abstract supertype for two-terminal voltage-source converter (VSC) HVDC formulations.
+Models per-terminal converters with quadratic / two-term losses
+(``a I^2 + b |I| + c``), a shared signed cable current, an explicit DC-side
+cable resistance (``v_f - v_t = (1/g) \\cdot I``), and (on AC networks) independent
+reactive-power control bounded by per-terminal PQ capability.
+"""
+abstract type AbstractTwoTerminalVSCFormulation <: AbstractTwoTerminalDCLineFormulation end
+
+"""
+Two-terminal VSC formulation that keeps the bilinear ``v \\cdot I`` and quadratic
+``I^2`` terms exact. Requires an NLP-capable solver (e.g. Ipopt).
+"""
+struct HVDCTwoTerminalVSCNLP <: AbstractTwoTerminalVSCFormulation end
+
+"""
+Two-terminal VSC formulation that uses SOS2 piecewise-linear surrogates for the
+bilinear ``v \\cdot I`` and quadratic ``I^2`` terms (so the loss model itself is
+mixed-integer linear) and enforces the per-terminal PQ capability via a linear
+outer-approximation of the disk ``p^2 + q^2 \\le \\text{rating}^2``: axis-aligned
+box constraints ``|p|, |q| \\le \\text{rating}`` always, plus four diagonal
+constraints ``|p| \\pm q \\le \\text{rating}\\sqrt{2}`` when the device-model
+attribute `use_octagon` (default `true`) is on. With the diagonals in place the
+feasible region is a regular octagon circumscribing the disk; turning them off
+leaves only the box.
+"""
+struct HVDCTwoTerminalVSCLP <: AbstractTwoTerminalVSCFormulation end
 
 ############################### AC/DC Converter Formulations #####################################
 abstract type AbstractConverterFormulation <: AbstractDeviceFormulation end
@@ -206,10 +231,10 @@ Abstract supertype for InterconnectingConverter formulations with quadratic loss
 abstract type AbstractQuadraticLossConverter <: AbstractConverterFormulation end
 
 """
-Quadratic Loss InterconnectingConverter using the Bin2 separable bilinear approximation
-(`v·i = ½((v+i)² − v² − i²)`) with a SOS2-based PWL approximation for x².
+Quadratic Loss InterconnectingConverter using the separable bilinear approximation
+(`v·i = ½((v+i)² − v² − i²)`) with a SOS2-based PWL approximation for x². Stays MILP.
 """
-struct Bin2QuadraticLossConverter <: AbstractQuadraticLossConverter end
+struct MILPQuadraticLossConverter <: AbstractQuadraticLossConverter end
 
 """
 Quadratic Loss InterconnectingConverter using exact bilinear (v·i) and quadratic (i²)
diff --git a/src/core/variables.jl b/src/core/variables.jl
index 8a41d0d..bb3ce26 100644
--- a/src/core/variables.jl
+++ b/src/core/variables.jl
@@ -157,22 +157,25 @@ Docs abbreviation: ``i_c^{dc}``
 struct ConverterCurrent <: VariableType end
 
 """
-Struct to dispatch the creation of DC Converter Positive Term Current Variables for DC formulations
-Docs abbreviation: ``i_c^{+,dc}``
+Positive part of an absolute-value current decomposition (``i = i^+ - i^-``).
+Used by both two-terminal HVDC links (cable current) and interconnecting converters.
+Docs abbreviation: ``i^{+,dc}``
 """
-struct ConverterPositiveCurrent <: VariableType end
+struct PositiveCurrent <: VariableType end
 
 """
-Struct to dispatch the creation of DC Converter Negative Term Current Variables for DC formulations
-Docs abbreviation: ``i_c^{-,dc}``
+Negative part of an absolute-value current decomposition (``i = i^+ - i^-``).
+Used by both two-terminal HVDC links (cable current) and interconnecting converters.
+Docs abbreviation: ``i^{-,dc}``
 """
-struct ConverterNegativeCurrent <: VariableType end
+struct NegativeCurrent <: VariableType end
 
 """
-Struct to dispatch the creation of DC Converter Binary for Absolute Value Current Variables for DC formulations
-Docs abbreviation: `\\nu_c``
+Binary indicator selecting which sign the current takes in the absolute-value
+decomposition (``i^+`` active when 1, ``i^-`` active when 0).
+Docs abbreviation: ``\\nu``
 """
-struct ConverterCurrentDirection <: VariableType end
+struct CurrentDirection <: VariableType end
 
 #########################################################
 #########################################################
@@ -384,6 +387,34 @@ Docs abbreviation: ``z``
 """
 struct HVDCPiecewiseBinaryLossVariable <: SparseVariableType end
 
+"""
+DC-side voltage at the from-terminal of a two-terminal HVDC link.
+Used by `HVDCTwoTerminalVSCNLP` formulations.
+Docs abbreviation: ``v_f^{dc}``
+"""
+struct HVDCFromDCVoltage <: VariableType end
+
+"""
+DC-side voltage at the to-terminal of a two-terminal HVDC link.
+Used by `HVDCTwoTerminalVSCNLP` formulations.
+Docs abbreviation: ``v_t^{dc}``
+"""
+struct HVDCToDCVoltage <: VariableType end
+
+"""
+Reactive power injected at the from-terminal AC bus by a two-terminal HVDC link.
+Added only when the formulation runs on an AC network model.
+Docs abbreviation: ``q_f``
+"""
+struct HVDCReactivePowerFromVariable <: VariableType end
+
+"""
+Reactive power injected at the to-terminal AC bus by a two-terminal HVDC link.
+Added only when the formulation runs on an AC network model.
+Docs abbreviation: ``q_t``
+"""
+struct HVDCReactivePowerToVariable <: VariableType end
+
 """
 Struct to dispatch the creation of Interface Flow Slack Up variables
 
diff --git a/src/mt_hvdc_models/HVDCsystems.jl b/src/mt_hvdc_models/HVDCsystems.jl
index b3fcde8..37dd418 100644
--- a/src/mt_hvdc_models/HVDCsystems.jl
+++ b/src/mt_hvdc_models/HVDCsystems.jl
@@ -120,38 +120,42 @@ end
 
 ## Binaries ###
 get_variable_binary(::Type{ConverterCurrent}, ::Type{PSY.InterconnectingConverter}, ::Type{<:AbstractConverterFormulation}) = false
-get_variable_binary(::Type{ConverterPositiveCurrent}, ::Type{PSY.InterconnectingConverter}, ::Type{<:AbstractConverterFormulation}) = false
-get_variable_binary(::Type{ConverterNegativeCurrent}, ::Type{PSY.InterconnectingConverter}, ::Type{<:AbstractConverterFormulation}) = false
-get_variable_binary(::Type{ConverterCurrentDirection}, ::Type{PSY.InterconnectingConverter}, ::Type{<:AbstractConverterFormulation}) = true
+get_variable_binary(::Type{PositiveCurrent}, ::Type{PSY.InterconnectingConverter}, ::Type{<:AbstractConverterFormulation}) = false
+get_variable_binary(::Type{NegativeCurrent}, ::Type{PSY.InterconnectingConverter}, ::Type{<:AbstractConverterFormulation}) = false
+get_variable_binary(::Type{CurrentDirection}, ::Type{PSY.InterconnectingConverter}, ::Type{<:AbstractConverterFormulation}) = true
 
 ### Warm Start ###
 get_variable_warm_start_value(::Type{ConverterCurrent}, d::PSY.InterconnectingConverter, ::Type{<:AbstractConverterFormulation}) = PSY.get_dc_current(d)
 
 ### Lower Bounds ###
 get_variable_lower_bound(::Type{ConverterCurrent}, d::PSY.InterconnectingConverter, ::Type{<:AbstractConverterFormulation}) = -PSY.get_max_dc_current(d)
-get_variable_lower_bound(::Type{ConverterPositiveCurrent}, d::PSY.InterconnectingConverter,::Type{<:AbstractConverterFormulation}) = 0.0
-get_variable_lower_bound(::Type{ConverterNegativeCurrent}, d::PSY.InterconnectingConverter,::Type{<:AbstractConverterFormulation}) = 0.0
+get_variable_lower_bound(::Type{PositiveCurrent}, d::PSY.InterconnectingConverter,::Type{<:AbstractConverterFormulation}) = 0.0
+get_variable_lower_bound(::Type{NegativeCurrent}, d::PSY.InterconnectingConverter,::Type{<:AbstractConverterFormulation}) = 0.0
 
 ### Upper Bounds ###
 get_variable_upper_bound(::Type{ConverterCurrent}, d::PSY.InterconnectingConverter, ::Type{<:AbstractConverterFormulation}) = PSY.get_max_dc_current(d)
-get_variable_upper_bound(::Type{ConverterPositiveCurrent}, d::PSY.InterconnectingConverter,::Type{<:AbstractConverterFormulation}) = PSY.get_max_dc_current(d)
-get_variable_upper_bound(::Type{ConverterNegativeCurrent}, d::PSY.InterconnectingConverter,::Type{<:AbstractConverterFormulation}) = PSY.get_max_dc_current(d)
+get_variable_upper_bound(::Type{PositiveCurrent}, d::PSY.InterconnectingConverter,::Type{<:AbstractConverterFormulation}) = PSY.get_max_dc_current(d)
+get_variable_upper_bound(::Type{NegativeCurrent}, d::PSY.InterconnectingConverter,::Type{<:AbstractConverterFormulation}) = PSY.get_max_dc_current(d)
 
 
+# Default `use_linear_loss = true`: the SOS2 PWL approximations already make the
+# model MIP, so the extra binary direction variable from the linear-loss
+# decomposition is free.
 function get_default_attributes(
     ::Type{PSY.InterconnectingConverter},
-    ::Type{Bin2QuadraticLossConverter},
+    ::Type{MILPQuadraticLossConverter},
 )
-    return Dict{String, Any}()
+    return Dict{String, Any}("use_linear_loss" => true)
 end
 
+# Default `use_linear_loss = false`: the formulation is a smooth NLP solvable by
+# Ipopt; enabling the linear-loss term introduces a binary direction variable
+# that pushes the model to MINLP, which pure NLP solvers cannot handle.
 function get_default_attributes(
     ::Type{PSY.InterconnectingConverter},
     ::Type{QuadraticLossConverter},
 )
-    return Dict{String, Any}(
-        "use_linear_loss" => false,
-    )
+    return Dict{String, Any}("use_linear_loss" => false)
 end
 
 #! format: on
@@ -540,20 +544,6 @@ end
 
 ############## Converters ##################
 
-_get_quadratic_term(loss_fn::PSY.QuadraticCurve) = PSY.get_quadratic_term(loss_fn)
-_get_quadratic_term(loss_fn) = 0.0
-
-_use_linear_loss(::Type{Bin2QuadraticLossConverter}, _) = true
-_use_linear_loss(::Type{QuadraticLossConverter}, model) =
-    get_attribute(model, "use_linear_loss")
-
-function _devices_with_linear_loss(devices)
-    return [
-        d for d in devices if
-        !iszero(PSY.get_proportional_term(PSY.get_loss_function(d)))
-    ]
-end
-
 function add_constraints!(
     container::OptimizationContainer,
     ::Type{ConverterLossConstraint},
@@ -570,10 +560,11 @@ function add_constraints!(
     vi_expr = get_expression(container, IOM.BilinearProductExpression, U, "vi")
     i_sq_expr = get_expression(container, IOM.QuadraticExpression, U, "i_sq")
     use_linear_loss =
-        _use_linear_loss(V, model) && !isempty(_devices_with_linear_loss(devices))
+        get_attribute(model, "use_linear_loss") &&
+        !isempty(_devices_with_linear_loss(devices))
     if use_linear_loss
-        i_pos_var = get_variable(container, ConverterPositiveCurrent, U)
-        i_neg_var = get_variable(container, ConverterNegativeCurrent, U)
+        i_pos_var = get_variable(container, PositiveCurrent, U)
+        i_neg_var = get_variable(container, NegativeCurrent, U)
     end
 
     ipc_names = [PSY.get_name(d) for d in devices]
@@ -589,66 +580,15 @@ function add_constraints!(
         b = PSY.get_proportional_term(loss_function)
         c = PSY.get_constant_term(loss_function)
         for t in time_steps
-            loss = a * i_sq_expr[name, t] + c
-            if use_linear_loss && !iszero(b)
-                loss += b * (i_pos_var[name, t] + i_neg_var[name, t])
-            end
-            loss_const[name, t] = JuMP.@constraint(
-                jump_model,
-                P_ac_var[name, t] == vi_expr[name, t] - loss
-            )
-        end
-    end
-    return
-end
-
-function add_constraints!(
-    container::OptimizationContainer,
-    ::Type{T},
-    devices::W,
-    ::DeviceModel{U, V},
-    ::NetworkModel{<:AbstractPowerModel},
-) where {
-    T <: CurrentAbsoluteValueConstraint,
-    U <: PSY.InterconnectingConverter,
-    V <: AbstractQuadraticLossConverter,
-    W <: Union{Vector{U}, IS.FlattenIteratorWrapper{U}},
-}
-    time_steps = get_time_steps(container)
-    names = [PSY.get_name(d) for d in devices]
-    jump_model = get_jump_model(container)
-    i_var = get_variable(container, ConverterCurrent, U)
-    i_pos_var = get_variable(container, ConverterPositiveCurrent, U)
-    i_neg_var = get_variable(container, ConverterNegativeCurrent, U)
-    i_dir_var = get_variable(container, ConverterCurrentDirection, U)
-
-    abs_val_const = add_constraints_container!(
-        container, CurrentAbsoluteValueConstraint, U, names, time_steps,
-    )
-    pos_ub_const = add_constraints_container!(
-        container, CurrentAbsoluteValueConstraint, U, names, time_steps;
-        meta = "pos_ub",
-    )
-    neg_ub_const = add_constraints_container!(
-        container, CurrentAbsoluteValueConstraint, U, names, time_steps;
-        meta = "neg_ub",
-    )
-
-    for d in devices
-        name = PSY.get_name(d)
-        i_max = PSY.get_max_dc_current(d)
-        for t in time_steps
-            abs_val_const[name, t] = JuMP.@constraint(
-                jump_model,
-                i_var[name, t] == i_pos_var[name, t] - i_neg_var[name, t]
-            )
-            pos_ub_const[name, t] = JuMP.@constraint(
-                jump_model,
-                i_pos_var[name, t] <= i_max * i_dir_var[name, t]
+            i_pos_t = use_linear_loss ? i_pos_var[name, t] : nothing
+            i_neg_t = use_linear_loss ? i_neg_var[name, t] : nothing
+            loss = _quadratic_converter_loss_expr(
+                a, b, c, i_sq_expr[name, t], i_pos_t, i_neg_t;
+                use_linear_loss = use_linear_loss,
             )
-            neg_ub_const[name, t] = JuMP.@constraint(
+            loss_const[name, t] = JuMP.@constraint(
                 jump_model,
-                i_neg_var[name, t] <= i_max * (1 - i_dir_var[name, t])
+                P_ac_var[name, t] == vi_expr[name, t] - loss,
             )
         end
     end
diff --git a/src/mt_hvdc_models/hvdcsystems_constructor.jl b/src/mt_hvdc_models/hvdcsystems_constructor.jl
index c1ec321..e8d00a3 100644
--- a/src/mt_hvdc_models/hvdcsystems_constructor.jl
+++ b/src/mt_hvdc_models/hvdcsystems_constructor.jl
@@ -54,15 +54,17 @@ function construct_device!(
     devices = get_available_components(model, sys)
     add_variables!(container, ActivePowerVariable, devices, T)
     add_variables!(container, ConverterCurrent, devices, T)
-    if _use_linear_loss(T, model)
-        ll_devices = _devices_with_linear_loss(devices)
-        if isempty(ll_devices)
-            @warn "use_linear_loss is enabled but every InterconnectingConverter has a zero proportional loss term; no linear-loss variables/constraints will be added."
-        else
-            add_variables!(container, ConverterPositiveCurrent, ll_devices, T)
-            add_variables!(container, ConverterNegativeCurrent, ll_devices, T)
-            add_variables!(container, ConverterCurrentDirection, ll_devices, T)
-        end
+
+    # use_linear_loss = true adds a binary direction variable (CurrentDirection).
+    # Both MILPQuadraticLossConverter and QuadraticLossConverter dispatch through
+    # here: on the NLP QuadraticLossConverter this pushes the model from a smooth
+    # NLP to MINLP (caller must supply a MINLP-capable solver); on
+    # MILPQuadraticLossConverter the SOS2 PWL approximations already make the
+    # model MILP, so the extra binary is free.
+    if get_attribute(model, "use_linear_loss")
+        _add_abs_value_decomposition_variables!(
+            container, devices, model, network_model,
+        )
     end
 
     add_to_expression!(
@@ -77,19 +79,19 @@ function construct_device!(
     return
 end
 
- function _voltage_expr_per_converter(
-     container::OptimizationContainer,
-     devices,
-     ipc_names::Vector{String},
-     time_steps,
- )
-     v_var = get_variable(container, DCVoltage, PSY.DCBus)
-     bus_names = [PSY.get_name(PSY.get_dc_bus(d)) for d in devices]
-     return JuMP.Containers.DenseAxisArray(
-         [v_var[b, t] for b in bus_names, t in time_steps],
-         ipc_names, time_steps,
-     )
- end
+function _voltage_expr_per_converter(
+    container::OptimizationContainer,
+    devices,
+    ipc_names::Vector{String},
+    time_steps,
+)
+    v_var = get_variable(container, DCVoltage, PSY.DCBus)
+    bus_names = [PSY.get_name(PSY.get_dc_bus(d)) for d in devices]
+    return JuMP.Containers.DenseAxisArray(
+        [v_var[b, t] for b in bus_names, t in time_steps],
+        ipc_names, time_steps,
+    )
+end
 
 function _converter_vi_bounds(devices)
     n = length(devices)
@@ -104,9 +106,11 @@ function _converter_vi_bounds(devices)
     return v_bounds, i_bounds
 end
 
-_quad_config(::Type{Bin2QuadraticLossConverter}) = IOM.SolverSOS2QuadConfig(DEFAULT_INTERPOLATION_LENGTH)
+_quad_config(::Type{MILPQuadraticLossConverter}) =
+    IOM.SolverSOS2QuadConfig(DEFAULT_INTERPOLATION_LENGTH)
 _quad_config(::Type{QuadraticLossConverter}) = IOM.NoQuadApproxConfig()
-_bilinear_config(::Type{Bin2QuadraticLossConverter}) = IOM.Bin2Config(IOM.SolverSOS2QuadConfig(DEFAULT_INTERPOLATION_LENGTH))
+_bilinear_config(::Type{MILPQuadraticLossConverter}) =
+    IOM.Bin2Config(IOM.SolverSOS2QuadConfig(DEFAULT_INTERPOLATION_LENGTH))
 _bilinear_config(::Type{QuadraticLossConverter}) = IOM.NoBilinearApproxConfig()
 
 function construct_device!(
@@ -148,20 +152,18 @@ function construct_device!(
         "vi",
     )
 
-    add_constraints!(container, ConverterLossConstraint, devices, model, network_model)
-    if _use_linear_loss(T, model)
-        ll_devices = _devices_with_linear_loss(devices)
-        if !isempty(ll_devices)
-            add_constraints!(
-                container,
-                CurrentAbsoluteValueConstraint,
-                ll_devices,
-                model,
-                network_model,
-            )
-        end
+    # PositiveCurrent / NegativeCurrent / CurrentDirection variables were added
+    # in the ArgumentConstructStage; only the decomposition constraints need
+    # the JuMP model now. ConverterLossConstraint reads these variables, so
+    # the decomposition constraints must be added first.
+    if get_attribute(model, "use_linear_loss")
+        _add_abs_value_decomposition_constraints!(
+            container, devices, model, network_model, ConverterCurrent,
+        )
     end
 
+    add_constraints!(container, ConverterLossConstraint, devices, model, network_model)
+
     add_feedforward_constraints!(container, model, devices)
     add_to_objective_function!(
         container, devices, model, get_network_formulation(network_model),
diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl
index 09f999e..46e8827 100644
--- a/src/network_models/pm_translator.jl
+++ b/src/network_models/pm_translator.jl
@@ -663,6 +663,27 @@ function get_branch_to_pm(
     return Dict{String, Any}()
 end
 
+function get_branch_to_pm(
+    ix::Int,
+    branch::PSY.TwoTerminalVSCLine,
+    ::Type{<:AbstractTwoTerminalVSCFormulation},
+    ::Type{<:AbstractPowerModel},
+)
+    return Dict{String, Any}()
+end
+
+# Trait: should this DeviceModel be skipped when translating two-terminal HVDC
+# branches into PowerModels? LCC and VSC formulations are handled by their own
+# constraints in POM rather than via PM's built-in DC line model, so they're
+# skipped here. Default is `false`; override per formulation via dispatch.
+_skip_pm_two_terminal_translation(::DeviceModel) = false
+_skip_pm_two_terminal_translation(
+    ::DeviceModel{PSY.TwoTerminalLCCLine, HVDCTwoTerminalLCC},
+) = true
+_skip_pm_two_terminal_translation(
+    ::DeviceModel{PSY.TwoTerminalVSCLine, <:AbstractTwoTerminalVSCFormulation},
+) = true
+
 function get_branches_to_pm(
     sys::PSY.System,
     network_model::NetworkModel{S},
@@ -720,10 +741,7 @@ function get_branches_to_pm(
     for (d, device_model) in branch_template
         comp_type = get_component_type(device_model)
         !(comp_type <: T) && continue
-        if comp_type <: PSY.TwoTerminalLCCLine &&
-           get_formulation(device_model) <: HVDCTwoTerminalLCC
-            continue
-        end
+        _skip_pm_two_terminal_translation(device_model) && continue
         start_idx += length(PM_branches)
         for (i, branch) in enumerate(get_available_components(device_model, sys))
             ix = i + start_idx
diff --git a/src/twoterminal_hvdc_models/TwoTerminalDC_branches.jl b/src/twoterminal_hvdc_models/TwoTerminalDC_branches.jl
index 5422b6f..d320236 100644
--- a/src/twoterminal_hvdc_models/TwoTerminalDC_branches.jl
+++ b/src/twoterminal_hvdc_models/TwoTerminalDC_branches.jl
@@ -1433,3 +1433,403 @@ function add_constraints!(
     end
     return
 end
+
+##############################################################################
+####################### Two-Terminal VSC Formulation #########################
+##############################################################################
+
+#! format: off
+
+# Variable trait methods for the shared cable current and DC voltages
+get_variable_binary(::Type{DCLineCurrentFlowVariable}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+get_variable_binary(::Type{HVDCFromDCVoltage}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+get_variable_binary(::Type{HVDCToDCVoltage}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+get_variable_binary(::Type{HVDCReactivePowerFromVariable}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+get_variable_binary(::Type{HVDCReactivePowerToVariable}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+get_variable_binary(::Type{PositiveCurrent}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+get_variable_binary(::Type{NegativeCurrent}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+get_variable_binary(::Type{CurrentDirection}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = true
+get_variable_binary(::Type{FlowActivePowerFromToVariable}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+get_variable_binary(::Type{FlowActivePowerToFromVariable}, ::Type{PSY.TwoTerminalVSCLine}, ::Type{<:AbstractTwoTerminalVSCFormulation}) = false
+
+# Warm starts
+get_variable_warm_start_value(::Type{DCLineCurrentFlowVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_dc_current(d)
+get_variable_warm_start_value(::Type{HVDCReactivePowerFromVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_reactive_power_from(d)
+get_variable_warm_start_value(::Type{HVDCReactivePowerToVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_reactive_power_to(d)
+get_variable_warm_start_value(::Type{FlowActivePowerFromToVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_active_power_flow(d)
+get_variable_warm_start_value(::Type{FlowActivePowerToFromVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = -PSY.get_active_power_flow(d)
+
+# Active power flow bounds (per-terminal)
+get_variable_lower_bound(::Type{FlowActivePowerFromToVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_active_power_limits_from(d).min
+get_variable_upper_bound(::Type{FlowActivePowerFromToVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_active_power_limits_from(d).max
+get_variable_lower_bound(::Type{FlowActivePowerToFromVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_active_power_limits_to(d).min
+get_variable_upper_bound(::Type{FlowActivePowerToFromVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_active_power_limits_to(d).max
+
+# Reactive power bounds (per-terminal)
+get_variable_lower_bound(::Type{HVDCReactivePowerFromVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_reactive_power_limits_from(d).min
+get_variable_upper_bound(::Type{HVDCReactivePowerFromVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_reactive_power_limits_from(d).max
+get_variable_lower_bound(::Type{HVDCReactivePowerToVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_reactive_power_limits_to(d).min
+get_variable_upper_bound(::Type{HVDCReactivePowerToVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_reactive_power_limits_to(d).max
+
+# DC voltage bounds (per-terminal)
+get_variable_lower_bound(::Type{HVDCFromDCVoltage}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_voltage_limits_from(d).min
+get_variable_upper_bound(::Type{HVDCFromDCVoltage}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_voltage_limits_from(d).max
+get_variable_lower_bound(::Type{HVDCToDCVoltage}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_voltage_limits_to(d).min
+get_variable_upper_bound(::Type{HVDCToDCVoltage}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = PSY.get_voltage_limits_to(d).max
+
+# Shared cable current bounds — must respect BOTH terminals' I_max ratings.
+# `_linear_loss_i_max(::PSY.TwoTerminalVSCLine)` lives in
+# common_models/quadratic_converter_loss.jl and returns the same `min(...)`.
+get_variable_lower_bound(::Type{DCLineCurrentFlowVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = -_linear_loss_i_max(d)
+get_variable_upper_bound(::Type{DCLineCurrentFlowVariable}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = _linear_loss_i_max(d)
+
+# Positive/negative parts: each in [0, i_max]
+get_variable_lower_bound(::Type{PositiveCurrent}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = 0.0
+get_variable_upper_bound(::Type{PositiveCurrent}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = _linear_loss_i_max(d)
+get_variable_lower_bound(::Type{NegativeCurrent}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = 0.0
+get_variable_upper_bound(::Type{NegativeCurrent}, d::PSY.TwoTerminalVSCLine, ::Type{<:AbstractTwoTerminalVSCFormulation}) = _linear_loss_i_max(d)
+
+#! format: on
+
+####################### VSC PQ-approx registration ###########################
+
+# Register the four IOM `QuadraticExpression` handles (`p_ft_sq`, `p_tf_sq`,
+# `q_f_sq`, `q_t_sq`) that the disk constraint reads. Only the NLP formulation
+# on an AC network actually emits the disk constraint, so only that combo
+# registers anything; the LP path and active-power-only networks no-op.
+function _register_pq_sq_expressions!(
+    container::OptimizationContainer,
+    devices,
+    line_names,
+    time_steps,
+    ::DeviceModel{PSY.TwoTerminalVSCLine, HVDCTwoTerminalVSCNLP},
+    ::NetworkModel{<:AbstractPowerModel},
+)
+    # This dispatch is the NLP path, so the quad config is fixed.
+    quad_cfg = IOM.NoQuadApproxConfig()
+    p_ft = get_variable(container, FlowActivePowerFromToVariable, PSY.TwoTerminalVSCLine)
+    p_tf = get_variable(container, FlowActivePowerToFromVariable, PSY.TwoTerminalVSCLine)
+    q_f = get_variable(container, HVDCReactivePowerFromVariable, PSY.TwoTerminalVSCLine)
+    q_t = get_variable(container, HVDCReactivePowerToVariable, PSY.TwoTerminalVSCLine)
+    p_ft_bounds = PSY.get_active_power_limits_from.(devices)
+    p_tf_bounds = PSY.get_active_power_limits_to.(devices)
+    q_f_bounds = PSY.get_reactive_power_limits_from.(devices)
+    q_t_bounds = PSY.get_reactive_power_limits_to.(devices)
+    IOM._add_quadratic_approx!(
+        quad_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps, p_ft, p_ft_bounds, "p_ft_sq",
+    )
+    IOM._add_quadratic_approx!(
+        quad_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps, p_tf, p_tf_bounds, "p_tf_sq",
+    )
+    IOM._add_quadratic_approx!(
+        quad_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps, q_f, q_f_bounds, "q_f_sq",
+    )
+    IOM._add_quadratic_approx!(
+        quad_cfg, container, PSY.TwoTerminalVSCLine,
+        line_names, time_steps, q_t, q_t_bounds, "q_t_sq",
+    )
+    return
+end
+
+# LP path: no disk constraint, so no p_sq/q_sq are needed.
+_register_pq_sq_expressions!(
+    ::OptimizationContainer, _devices, _names, _times,
+    ::DeviceModel{PSY.TwoTerminalVSCLine, HVDCTwoTerminalVSCLP},
+    ::NetworkModel,
+) = nothing
+
+# Active-power-only networks don't carry reactive variables at all.
+_register_pq_sq_expressions!(
+    ::OptimizationContainer, _devices, _names, _times,
+    ::DeviceModel{PSY.TwoTerminalVSCLine, HVDCTwoTerminalVSCNLP},
+    ::NetworkModel{<:AbstractActivePowerModel},
+) = nothing
+
+####################### VSC core constraints ################################
+
+# Cable Ohm's law:  v_f - v_t = (1/g) * I
+function add_constraints!(
+    container::OptimizationContainer,
+    ::Type{HVDCCableOhmsLawConstraint},
+    devices::Union{Vector{U}, IS.FlattenIteratorWrapper{U}},
+    ::DeviceModel{U, F},
+    ::NetworkModel{<:AbstractPowerModel},
+) where {U <: PSY.TwoTerminalVSCLine, F <: AbstractTwoTerminalVSCFormulation}
+    time_steps = get_time_steps(container)
+    names = [PSY.get_name(d) for d in devices]
+    jump_model = get_jump_model(container)
+    v_f = get_variable(container, HVDCFromDCVoltage, U)
+    v_t = get_variable(container, HVDCToDCVoltage, U)
+    i_var = get_variable(container, DCLineCurrentFlowVariable, U)
+
+    cons = add_constraints_container!(
+        container, HVDCCableOhmsLawConstraint, U, names, time_steps,
+    )
+
+    for d in devices
+        name = PSY.get_name(d)
+        g = PSY.get_g(d)
+        for t in time_steps
+            cons[name, t] = if iszero(g)
+                JuMP.@constraint(jump_model, i_var[name, t] == 0)
+            else
+                JuMP.@constraint(
+                    jump_model,
+                    v_f[name, t] - v_t[name, t] == (1.0 / g) * i_var[name, t],
+                )
+            end
+        end
+    end
+    return
+end
+
+# Per-terminal converter power balance:
+#   p_ft ==  v_f * I + (a_f * I^2 + b_f * |I| + c_f)
+#   p_tf == -v_t * I + (a_t * I^2 + b_t * |I| + c_t)
+# Active power enters `ActivePowerBalance` with a -1 multiplier (the AC bus
+# sees the converter as a load drawing p_ft / p_tf), so positive values of
+# `FlowActivePower*Variable` correspond to power flowing AC → DC at the
+# respective terminal.
+function add_constraints!(
+    container::OptimizationContainer,
+    ::Type{HVDCVSCConverterPowerConstraint},
+    devices::Union{Vector{U}, IS.FlattenIteratorWrapper{U}},
+    model::DeviceModel{U, F},
+    ::NetworkModel{<:AbstractPowerModel},
+) where {U <: PSY.TwoTerminalVSCLine, F <: AbstractTwoTerminalVSCFormulation}
+    time_steps = get_time_steps(container)
+    names = [PSY.get_name(d) for d in devices]
+    jump_model = get_jump_model(container)
+
+    p_ft = get_variable(container, FlowActivePowerFromToVariable, U)
+    p_tf = get_variable(container, FlowActivePowerToFromVariable, U)
+    vi_expr_ft = get_expression(container, IOM.BilinearProductExpression, U, "vi_ft")
+    vi_expr_tf = get_expression(container, IOM.BilinearProductExpression, U, "vi_tf")
+    i_sq_expr = get_expression(container, IOM.QuadraticExpression, U, "i_sq")
+
+    use_linear_loss =
+        get_attribute(model, "use_linear_loss") &&
+        !isempty(_devices_with_linear_loss(devices))
+    if use_linear_loss
+        i_pos_var = get_variable(container, PositiveCurrent, U)
+        i_neg_var = get_variable(container, NegativeCurrent, U)
+    end
+
+    cons_ft = add_constraints_container!(
+        container, HVDCVSCConverterPowerConstraint, U, names, time_steps; meta = "ft",
+    )
+    cons_tf = add_constraints_container!(
+        container, HVDCVSCConverterPowerConstraint, U, names, time_steps; meta = "tf",
+    )
+
+    for d in devices
+        name = PSY.get_name(d)
+        loss_from = PSY.get_converter_loss_from(d)
+        loss_to = PSY.get_converter_loss_to(d)
+        a_f = _get_quadratic_term(loss_from)
+        b_f = PSY.get_proportional_term(loss_from)
+        c_f = PSY.get_constant_term(loss_from)
+        a_t = _get_quadratic_term(loss_to)
+        b_t = PSY.get_proportional_term(loss_to)
+        c_t = PSY.get_constant_term(loss_to)
+        for t in time_steps
+            i_pos_t = use_linear_loss ? i_pos_var[name, t] : nothing
+            i_neg_t = use_linear_loss ? i_neg_var[name, t] : nothing
+            loss_ft = _quadratic_converter_loss_expr(
+                a_f, b_f, c_f, i_sq_expr[name, t], i_pos_t, i_neg_t;
+                use_linear_loss = use_linear_loss,
+            )
+            loss_tf = _quadratic_converter_loss_expr(
+                a_t, b_t, c_t, i_sq_expr[name, t], i_pos_t, i_neg_t;
+                use_linear_loss = use_linear_loss,
+            )
+            cons_ft[name, t] = JuMP.@constraint(
+                jump_model,
+                p_ft[name, t] == vi_expr_ft[name, t] + loss_ft,
+            )
+            cons_tf[name, t] = JuMP.@constraint(
+                jump_model,
+                p_tf[name, t] == -vi_expr_tf[name, t] + loss_tf,
+            )
+        end
+    end
+    return
+end
+
+# PQ capability — exact disk for the NLP formulation. `p_*_sq` / `q_*_sq` are
+# the `IOM.QuadraticExpression` handles registered by
+# `_register_pq_sq_expressions!`. Under `NoQuadApproxConfig` (what
+# `HVDCTwoTerminalVSCNLP` uses) they are exact QuadExprs, so the constraint is
+# the smooth `p² + q² ≤ s²` and the model stays an NLP.
+function add_constraints!(
+    container::OptimizationContainer,
+    ::Type{HVDCVSCApparentPowerLimitConstraint},
+    devices::Union{Vector{U}, IS.FlattenIteratorWrapper{U}},
+    ::DeviceModel{U, HVDCTwoTerminalVSCNLP},
+    ::NetworkModel{<:AbstractPowerModel},
+) where {U <: PSY.TwoTerminalVSCLine}
+    time_steps = get_time_steps(container)
+    names = [PSY.get_name(d) for d in devices]
+    jump_model = get_jump_model(container)
+
+    p_ft_sq = get_expression(container, IOM.QuadraticExpression, U, "p_ft_sq")
+    p_tf_sq = get_expression(container, IOM.QuadraticExpression, U, "p_tf_sq")
+    q_f_sq = get_expression(container, IOM.QuadraticExpression, U, "q_f_sq")
+    q_t_sq = get_expression(container, IOM.QuadraticExpression, U, "q_t_sq")
+
+    cons_f = add_constraints_container!(
+        container, HVDCVSCApparentPowerLimitConstraint, U, names, time_steps;
+        meta = "from",
+    )
+    cons_t = add_constraints_container!(
+        container, HVDCVSCApparentPowerLimitConstraint, U, names, time_steps;
+        meta = "to",
+    )
+
+    for d in devices
+        name = PSY.get_name(d)
+        s_f2 = PSY.get_rating_from(d)^2
+        s_t2 = PSY.get_rating_to(d)^2
+        for t in time_steps
+            cons_f[name, t] = JuMP.@constraint(
+                jump_model, p_ft_sq[name, t] + q_f_sq[name, t] <= s_f2,
+            )
+            cons_t[name, t] = JuMP.@constraint(
+                jump_model, p_tf_sq[name, t] + q_t_sq[name, t] <= s_t2,
+            )
+        end
+    end
+    return
+end
+
+# PQ capability — linear outer-approximation for the LP formulation.
+#
+# We always add the axis-aligned box  |p|, |q| ≤ rating.  When the
+# device-model attribute `use_octagon` (default `true`) is on, we also add
+# the four 45°-rotated diagonals  |p| ± q ≤ rating·√2 ; their intersection
+# with the box is a regular octagon circumscribing the disk
+# p² + q² ≤ rating².
+#
+# Outer-approximation proof: for any (p, q) on the disk, p² ≤ p² + q² ≤ r²
+# gives |p|, |q| ≤ r, and Cauchy–Schwarz gives (|p|+|q|)² ≤ 2(p²+q²) ≤ 2r²
+# so |p|+|q| ≤ r√2. Both half-plane families contain the disk, and so does
+# their intersection. The octagon is loose by at most ≈8.2% in area
+# (octagon-to-disk area ratio 8·tan(π/8)/π ≈ 1.082).
+function add_constraints!(
+    container::OptimizationContainer,
+    ::Type{HVDCVSCApparentPowerLimitConstraint},
+    devices::Union{Vector{U}, IS.FlattenIteratorWrapper{U}},
+    model::DeviceModel{U, HVDCTwoTerminalVSCLP},
+    ::NetworkModel{<:AbstractPowerModel},
+) where {U <: PSY.TwoTerminalVSCLine}
+    time_steps = get_time_steps(container)
+    names = [PSY.get_name(d) for d in devices]
+    jump_model = get_jump_model(container)
+
+    p_ft = get_variable(container, FlowActivePowerFromToVariable, U)
+    p_tf = get_variable(container, FlowActivePowerToFromVariable, U)
+    q_f = get_variable(container, HVDCReactivePowerFromVariable, U)
+    q_t = get_variable(container, HVDCReactivePowerToVariable, U)
+
+    use_octagon = get_attribute(model, "use_octagon")
+    side_tags = if use_octagon
+        ("from_p_ub", "from_p_lb", "from_q_ub", "from_q_lb",
+            "to_p_ub", "to_p_lb", "to_q_ub", "to_q_lb",
+            "from_pp", "from_pn", "from_np", "from_nn",
+            "to_pp", "to_pn", "to_np", "to_nn")
+    else
+        ("from_p_ub", "from_p_lb", "from_q_ub", "from_q_lb",
+            "to_p_ub", "to_p_lb", "to_q_ub", "to_q_lb")
+    end
+    cons = Dict{String, Any}()
+    for tag in side_tags
+        cons[tag] = add_constraints_container!(
+            container, HVDCVSCApparentPowerLimitConstraint, U,
+            names, time_steps; meta = tag,
+        )
+    end
+
+    side_specs = (
+        (prefix = "from", p_var = p_ft, q_var = q_f, rating_getter = PSY.get_rating_from),
+        (prefix = "to", p_var = p_tf, q_var = q_t, rating_getter = PSY.get_rating_to),
+    )
+    for d in devices
+        name = PSY.get_name(d)
+        for spec in side_specs
+            rating = spec.rating_getter(d)
+            diag = rating * sqrt(2.0)
+            prefix = spec.prefix
+            p_var, q_var = spec.p_var, spec.q_var
+            for t in time_steps
+                cons[prefix * "_p_ub"][name, t] =
+                    JuMP.@constraint(jump_model, p_var[name, t] <= rating)
+                cons[prefix * "_p_lb"][name, t] =
+                    JuMP.@constraint(jump_model, -p_var[name, t] <= rating)
+                cons[prefix * "_q_ub"][name, t] =
+                    JuMP.@constraint(jump_model, q_var[name, t] <= rating)
+                cons[prefix * "_q_lb"][name, t] =
+                    JuMP.@constraint(jump_model, -q_var[name, t] <= rating)
+                if use_octagon
+                    cons[prefix * "_pp"][name, t] =
+                        JuMP.@constraint(
+                            jump_model,
+                            p_var[name, t] + q_var[name, t] <= diag
+                        )
+                    cons[prefix * "_pn"][name, t] =
+                        JuMP.@constraint(
+                            jump_model,
+                            p_var[name, t] - q_var[name, t] <= diag
+                        )
+                    cons[prefix * "_np"][name, t] =
+                        JuMP.@constraint(
+                            jump_model,
+                            -p_var[name, t] + q_var[name, t] <= diag
+                        )
+                    cons[prefix * "_nn"][name, t] =
+                        JuMP.@constraint(
+                            jump_model,
+                            -p_var[name, t] - q_var[name, t] <= diag
+                        )
+                end
+            end
+        end
+    end
+    return
+end
+
+####################### VSC defaults #########################################
+
+function get_default_time_series_names(
+    ::Type{PSY.TwoTerminalVSCLine},
+    ::Type{<:AbstractTwoTerminalVSCFormulation},
+)
+    return Dict{Type{<:TimeSeriesParameter}, String}()
+end
+
+# Default `use_linear_loss = false`: HVDCTwoTerminalVSCNLP is a smooth NLP solvable
+# by Ipopt; enabling the linear-loss term introduces a binary direction variable
+# that pushes the model to MINLP, which pure NLP solvers cannot handle.
+function get_default_attributes(
+    ::Type{PSY.TwoTerminalVSCLine},
+    ::Type{HVDCTwoTerminalVSCNLP},
+)
+    return Dict{String, Any}("use_linear_loss" => false)
+end
+
+# `use_linear_loss = true`: the SOS2 PWL loss model already pushes the
+# formulation into MILP territory, so the extra binary direction variable from
+# the linear-loss decomposition is free.
+# `use_octagon = true`: adds the four diagonals |p| ± q ≤ rating·√2 on top of
+# the axis-aligned box |p|, |q| ≤ rating. The intersection is a regular octagon
+# circumscribing the disk p² + q² ≤ rating² and is a guaranteed outer
+# approximation (loose by at most ≈8.2% in area). Setting it to false leaves
+# only the box, which is cheaper but a looser linear envelope of the disk.
+function get_default_attributes(
+    ::Type{PSY.TwoTerminalVSCLine},
+    ::Type{HVDCTwoTerminalVSCLP},
+)
+    return Dict{String, Any}("use_linear_loss" => true, "use_octagon" => true)
+end
diff --git a/test/Project.toml b/test/Project.toml
index a7425a5..43494b2 100644
--- a/test/Project.toml
+++ b/test/Project.toml
@@ -35,7 +35,7 @@ UUIDs = "cf7118a7-6976-5b1a-9a39-7adc72f591a4"
 InfrastructureSystems = {rev = "IS4", url = "https://github.com/Sienna-Platform/InfrastructureSystems.jl"}
 PowerSystemCaseBuilder = {rev = "psy6", url = "https://github.com/Sienna-Platform/PowerSystemCaseBuilder.jl"}
 PowerSystems = {rev = "psy6", url = "https://github.com/Sienna-Platform/PowerSystems.jl"}
-
+InfrastructureOptimizationModels = {rev = "ac/hvdc-vsc", url = "https://github.com/Sienna-Platform/InfrastructureOptimizationModels.jl"}
 
 [compat]
 HiGHS = "1"
diff --git a/test/runtests.jl b/test/runtests.jl
index 7d8fc42..d290015 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -22,8 +22,14 @@ const DISABLED_TEST_FILES = [  # Can generate with ls -1 test | grep "test_.*.jl
 # "test_device_synchronous_condenser_constructors.jl",
 # "test_device_thermal_generation_constructors.jl",
 # "test_formulation_combinations.jl",
+# "test_import_export_cost.jl",
 # "test_initialization_problem.jl",
+# "test_is_time_variant_proportional.jl",
+# "test_market_bid_cost.jl",
+# "test_mbc_parameter_population.jl",
 # "test_model_decision.jl",
+# "test_multi_interval.jl",
+# "test_network_constructors.jl",
 # "test_network_constructors_with_dlr.jl",
 # "test_problem_template.jl",
 # "test_storage_device_models.jl",
diff --git a/test/test_device_hvdc.jl b/test/test_device_hvdc.jl
index 62652e8..2632e26 100644
--- a/test/test_device_hvdc.jl
+++ b/test/test_device_hvdc.jl
@@ -97,7 +97,7 @@ end
     @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
 end
 
-@testset "HVDC System with Losses Network (Bin2QuadraticLossConverter)" begin
+@testset "HVDC System with Losses Network (MILPQuadraticLossConverter)" begin
     sys = _generate_test_hvdc_sys()
     template = OperationsProblemTemplate()
     set_device_model!(template, ThermalStandard, ThermalDispatchNoMin)
@@ -106,7 +106,7 @@ end
     set_device_model!(template, TModelHVDCLine, DCLossyLine)
     ipc_model = DeviceModel(
         InterconnectingConverter,
-        Bin2QuadraticLossConverter,
+        MILPQuadraticLossConverter,
     )
     set_device_model!(template, ipc_model)
     set_hvdc_network_model!(template, VoltageDispatchHVDCNetworkModel)
@@ -147,7 +147,7 @@ end
     @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
 end
 
-@testset "HVDC Bin2 vs NLP QuadraticLossConverter agreement" begin
+@testset "HVDC MILP vs NLP QuadraticLossConverter agreement" begin
     function _build_and_solve(formulation, optimizer)
         sys = _generate_test_hvdc_sys()
         template = OperationsProblemTemplate()
@@ -176,41 +176,187 @@ end
         return model
     end
 
-    bin2_model = _build_and_solve(Bin2QuadraticLossConverter, HiGHS_optimizer)
+    milp_model = _build_and_solve(MILPQuadraticLossConverter, HiGHS_optimizer)
     nlp_model = _build_and_solve(QuadraticLossConverter, ipopt_optimizer)
 
-    bin2_obj = IOM.get_objective_value(OptimizationProblemOutputs(bin2_model))
+    milp_obj = IOM.get_objective_value(OptimizationProblemOutputs(milp_model))
     nlp_obj = IOM.get_objective_value(OptimizationProblemOutputs(nlp_model))
 
-    # Bin2 is a relaxation/PWL approximation of the exact NLP. The two objectives
-    # should agree to within a few percent on this small system.
-    @test isapprox(bin2_obj, nlp_obj; rtol = 0.05)
+    @test isapprox(milp_obj, nlp_obj; rtol = 0.05)
 end
 
-@testset "HVDC linear-loss warning when all converters have b=0" begin
-    sys = _generate_test_hvdc_sys()
-    for ipc in get_components(InterconnectingConverter, sys)
-        set_loss_function!(ipc, QuadraticCurve(0.01, 0.0, 0.0))
-    end
-    template = OperationsProblemTemplate()
+##############################################################################
+################ Two-Terminal VSC HVDC tests #################################
+##############################################################################
+
+# Build a small AC test system and replace an AC line with a TwoTerminalVSCLine
+# so we have a concrete VSC device to exercise the formulation against.
+function _generate_test_vsc_sys(;
+    g = 50.0,
+    rating_from = 2.0,
+    rating_to = 2.0,
+    loss_a = 0.01,
+    loss_b = 0.0,
+    loss_c = 0.0,
+)
+    sys = build_system(PSITestSystems, "c_sys5_uc"; force_build = true)
+    line = get_component(Line, sys, "1")
+    remove_component!(sys, line)
+
+    vsc = TwoTerminalVSCLine(;
+        name = get_name(line),
+        available = true,
+        arc = get_arc(line),
+        active_power_flow = 0.0,
+        rating = max(rating_from, rating_to),
+        active_power_limits_from = (min = -rating_from, max = rating_from),
+        active_power_limits_to = (min = -rating_to, max = rating_to),
+        g = g,
+        dc_current = 0.0,
+        reactive_power_from = 0.0,
+        dc_voltage_control_from = true,
+        ac_voltage_control_from = true,
+        dc_setpoint_from = 0.0,
+        ac_setpoint_from = 1.0,
+        converter_loss_from = QuadraticCurve(loss_a, loss_b, loss_c),
+        max_dc_current_from = 5.0,
+        rating_from = rating_from,
+        reactive_power_limits_from = (min = -rating_from, max = rating_from),
+        power_factor_weighting_fraction_from = 1.0,
+        voltage_limits_from = (min = 0.95, max = 1.05),
+        reactive_power_to = 0.0,
+        dc_voltage_control_to = true,
+        ac_voltage_control_to = true,
+        dc_setpoint_to = 0.0,
+        ac_setpoint_to = 1.0,
+        converter_loss_to = QuadraticCurve(loss_a, loss_b, loss_c),
+        max_dc_current_to = 5.0,
+        rating_to = rating_to,
+        reactive_power_limits_to = (min = -rating_to, max = rating_to),
+        power_factor_weighting_fraction_to = 1.0,
+        voltage_limits_to = (min = 0.95, max = 1.05),
+    )
+    add_component!(sys, vsc)
+    return sys
+end
+
+function _build_vsc_model(formulation, network, optimizer; sys = _generate_test_vsc_sys())
+    template = OperationsProblemTemplate(NetworkModel(network))
     set_device_model!(template, ThermalStandard, ThermalDispatchNoMin)
+    set_device_model!(template, RenewableDispatch, RenewableFullDispatch)
     set_device_model!(template, PowerLoad, StaticPowerLoad)
     set_device_model!(template, DeviceModel(Line, StaticBranch))
-    set_device_model!(template, TModelHVDCLine, DCLossyLine)
-    set_device_model!(
-        template,
-        DeviceModel(InterconnectingConverter, Bin2QuadraticLossConverter),
-    )
-    set_hvdc_network_model!(template, VoltageDispatchHVDCNetworkModel)
-    model = DecisionModel(
-        template,
-        sys;
-        store_variable_names = true,
-        optimizer = HiGHS_optimizer,
-    )
-    @test_logs(
-        (:warn, r"every InterconnectingConverter has a zero proportional loss"),
-        match_mode = :any,
-        build!(model; output_dir = mktempdir(; cleanup = true)),
+    set_device_model!(template, DeviceModel(TwoTerminalVSCLine, formulation))
+    return DecisionModel(
+        template, sys; store_variable_names = true, optimizer = optimizer,
     )
 end
+
+@testset "HVDC Two-Terminal VSC (NLP) on DCP" begin
+    model = _build_vsc_model(HVDCTwoTerminalVSCNLP, DCPPowerModel, ipopt_optimizer)
+    @test build!(model; output_dir = mktempdir(; cleanup = true)) ==
+          IOM.ModelBuildStatus.BUILT
+    @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
+end
+
+@testset "HVDC Two-Terminal VSC on AC (NLP)" begin
+    model = _build_vsc_model(HVDCTwoTerminalVSCNLP, ACPPowerModel, ipopt_optimizer)
+    @test build!(model; output_dir = mktempdir(; cleanup = true)) ==
+          IOM.ModelBuildStatus.BUILT
+    @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
+end
+
+@testset "HVDC Two-Terminal VSC (LP) on DCP" begin
+    model = _build_vsc_model(HVDCTwoTerminalVSCLP, DCPPowerModel, HiGHS_optimizer)
+    @test build!(model; output_dir = mktempdir(; cleanup = true)) ==
+          IOM.ModelBuildStatus.BUILT
+    @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
+end
+
+# Omitted: HVDCTwoTerminalVSCLP on ACPPowerModel — HiGHS cannot solve the
+# ACP network's trigonometric (cos/sin) branch ohms-law constraints; we'd
+# need an MINLP-capable solver (e.g. Gurobi) for that combination.
+
+@testset "HVDC VSC LP vs NLP objective agreement" begin
+    # On a DC network the PQ disk constraint is inactive (no reactive
+    # variables exist), so the LP and NLP differ only by the i² loss model
+    # (SOS2 PWL vs exact). For a smooth convex loss curve the two should agree
+    # within a few percent.
+    function _solve(formulation, optimizer)
+        sys = _generate_test_vsc_sys()
+        model = _build_vsc_model(formulation, DCPPowerModel, optimizer; sys = sys)
+        @test build!(model; output_dir = mktempdir(; cleanup = true)) ==
+              IOM.ModelBuildStatus.BUILT
+        @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
+        return IOM.get_optimization_container(model).optimizer_stats.objective_value
+    end
+    lp_obj = _solve(HVDCTwoTerminalVSCLP, HiGHS_optimizer)
+    nlp_obj = _solve(HVDCTwoTerminalVSCNLP, ipopt_optimizer)
+    @test isapprox(lp_obj, nlp_obj; rtol = 0.05)
+end
+
+@testset "HVDC VSC LP: octagon is at least as tight as box-only" begin
+    # With `use_octagon = true` the LP adds four diagonals on top of the box,
+    # shrinking the feasible region (or leaving it unchanged if the diagonals
+    # are non-binding). The min-cost objective should therefore be ≥ the
+    # box-only case.
+    function _solve_lp(use_octagon)
+        sys = _generate_test_vsc_sys()
+        template = OperationsProblemTemplate(NetworkModel(DCPPowerModel))
+        set_device_model!(template, ThermalStandard, ThermalDispatchNoMin)
+        set_device_model!(template, RenewableDispatch, RenewableFullDispatch)
+        set_device_model!(template, PowerLoad, StaticPowerLoad)
+        set_device_model!(template, DeviceModel(Line, StaticBranch))
+        set_device_model!(
+            template,
+            DeviceModel(
+                TwoTerminalVSCLine,
+                HVDCTwoTerminalVSCLP;
+                attributes = Dict("use_linear_loss" => true, "use_octagon" => use_octagon),
+            ),
+        )
+        model = DecisionModel(
+            template, sys; store_variable_names = true, optimizer = HiGHS_optimizer,
+        )
+        @test build!(model; output_dir = mktempdir(; cleanup = true)) ==
+              IOM.ModelBuildStatus.BUILT
+        @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
+        return IOM.get_optimization_container(model).optimizer_stats.objective_value
+    end
+    box_only_obj = _solve_lp(false)
+    octagon_obj = _solve_lp(true)
+    @test octagon_obj >= box_only_obj - 1e-6
+end
+
+@testset "HVDC VSC: higher cable resistance increases cost" begin
+    # Smaller g => larger R = 1/g => more losses => optimum should not improve.
+    function _solve_with_g(g_value)
+        sys = _generate_test_vsc_sys(; g = g_value)
+        model = _build_vsc_model(
+            HVDCTwoTerminalVSCLP, DCPPowerModel, HiGHS_optimizer; sys = sys,
+        )
+        @test build!(model; output_dir = mktempdir(; cleanup = true)) ==
+              IOM.ModelBuildStatus.BUILT
+        @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
+        return IOM.get_optimization_container(model).optimizer_stats.objective_value
+    end
+    low_R_obj = _solve_with_g(100.0)   # large g, small R
+    high_R_obj = _solve_with_g(20.0)   # smaller g, larger R
+    @test high_R_obj >= low_R_obj - 1e-6
+end
+
+@testset "HVDC VSC: tighter PQ rating raises cost on AC" begin
+    function _solve_with_rating(s)
+        sys = _generate_test_vsc_sys(; rating_from = s, rating_to = s)
+        model = _build_vsc_model(
+            HVDCTwoTerminalVSCNLP, ACPPowerModel, ipopt_optimizer; sys = sys,
+        )
+        @test build!(model; output_dir = mktempdir(; cleanup = true)) ==
+              IOM.ModelBuildStatus.BUILT
+        @test solve!(model) == IOM.RunStatus.SUCCESSFULLY_FINALIZED
+        return IOM.get_optimization_container(model).optimizer_stats.objective_value
+    end
+    looser = _solve_with_rating(2.0)
+    tighter = _solve_with_rating(1.0)
+    @test tighter >= looser - 1e-6
+end