diff --git a/src/Ybus.jl b/src/Ybus.jl index fe17d043d..3c8dc5f20 100644 --- a/src/Ybus.jl +++ b/src/Ybus.jl @@ -435,9 +435,104 @@ end _get_tap(::PSY.Transformer2W) = one(YBUS_ELTYPE) _get_tap(br::PSY.TwoWindingTransformer) = PSY.get_tap(br) +function _interpolate_correction_factor(curve::IS.PiecewiseLinearData, x::Real) + points = IS.get_points(curve) + x = clamp(x, points[1].x, points[end].x) + for i in 1:(length(points) - 1) + if x <= points[i + 1].x + dx = points[i + 1].x - points[i].x + iszero(dx) && return points[i].y + t = (x - points[i].x) / dx + return points[i].y + t * (points[i + 1].y - points[i].y) + end + end + return points[end].y +end + +function _compute_icd_factor( + br::PSY.TwoWindingTransformer, + ict::PSY.ImpedanceCorrectionData, +) + curve = PSY.get_impedance_correction_curve(ict) + mode = PSY.get_transformer_control_mode(ict) + x = if mode == PSY.ImpedanceCorrectionTransformerControlMode.TAP_RATIO + abs(_get_tap(br)) + else # PHASE_SHIFT_ANGLE — table x-values are in degrees; α is stored in radians + rad2deg(PSY.get_α(br)) + end + return _interpolate_correction_factor(curve, x) +end + +function _compute_icd_factor( + br::PSY.ThreeWindingTransformer, + ict::PSY.ImpedanceCorrectionData, + winding_num::Int, +) + curve = PSY.get_impedance_correction_curve(ict) + mode = PSY.get_transformer_control_mode(ict) + x = if mode == PSY.ImpedanceCorrectionTransformerControlMode.TAP_RATIO + abs( + ( + PSY.get_primary_turns_ratio(br), + PSY.get_secondary_turns_ratio(br), + PSY.get_tertiary_turns_ratio(br), + )[winding_num], + ) + else # PHASE_SHIFT_ANGLE — table x-values are in degrees; angles stored in radians + rad2deg( + ( + PSY.get_α_primary(br), + PSY.get_α_secondary(br), + PSY.get_α_tertiary(br), + )[winding_num], + ) + end + return _interpolate_correction_factor(curve, x) +end + +function _winding_category_to_number(cat::PSY.WindingCategoryModule.WindingCategory)::Int + cat == PSY.WindingCategory.PRIMARY_WINDING && return 1 + cat == PSY.WindingCategory.SECONDARY_WINDING && return 2 + return 3 +end + +function _build_icd_correction_map(sys::PSY.System)::Union{Nothing, ICDCorrectionMap} + pairs_2w = PSY.get_component_supplemental_attribute_pairs( + PSY.TwoWindingTransformer, + PSY.ImpedanceCorrectionData, + sys, + ) + pairs_3w = PSY.get_component_supplemental_attribute_pairs( + PSY.ThreeWindingTransformer, + PSY.ImpedanceCorrectionData, + sys, + ) + isempty(pairs_2w) && isempty(pairs_3w) && return nothing + map_2w = Dict{Base.UUID, Float64}() + for pair in pairs_2w + map_2w[IS.get_uuid(pair.component)] = + _compute_icd_factor(pair.component, pair.supplemental_attribute) + end + map_3w = Dict{Tuple{Base.UUID, Int}, Float64}() + for pair in pairs_3w + wnum = _winding_category_to_number( + PSY.get_transformer_winding(pair.supplemental_attribute), + ) + map_3w[(IS.get_uuid(pair.component), wnum)] = + _compute_icd_factor(pair.component, pair.supplemental_attribute, wnum) + end + return ICDCorrectionMap(map_2w, map_3w) +end + +function _get_impedance_correction_factor(br::PSY.TwoWindingTransformer) + attrs = PSY.get_supplemental_attributes(PSY.ImpedanceCorrectionData, br) + isempty(attrs) && return one(Float64) + return _compute_icd_factor(br, only(attrs)) +end + """Ybus entries for a `Transformer2W`, `TapTransformer`, or `PhaseShiftingTransformer`.""" -function ybus_branch_entries(br::PSY.TwoWindingTransformer) - Y_t = 1 / (PSY.get_r(br) + PSY.get_x(br) * 1im) +function ybus_branch_entries(br::PSY.TwoWindingTransformer, correction::Float64) + Y_t = 1 / ((PSY.get_r(br) + PSY.get_x(br) * 1im) * correction) tap = _get_tap(br) * exp(PSY.get_α(br) * 1im) c_tap = _get_tap(br) * exp(-1 * PSY.get_α(br) * 1im) y_shunt = PSY.get_primary_shunt(br) @@ -453,12 +548,32 @@ function ybus_branch_entries(br::PSY.TwoWindingTransformer) return (Y11 + y_shunt, Y12, Y21, Y22) end +function ybus_branch_entries(br::PSY.TwoWindingTransformer) + return ybus_branch_entries(br, _get_impedance_correction_factor(br)) +end + +function _get_impedance_correction_factor(tp::ThreeWindingTransformerWinding) + br = get_transformer(tp) + winding_num = get_winding_number(tp) + target = ( + PSY.WindingCategory.PRIMARY_WINDING, + PSY.WindingCategory.SECONDARY_WINDING, + PSY.WindingCategory.TERTIARY_WINDING, + )[winding_num] + attrs = PSY.get_supplemental_attributes(PSY.ImpedanceCorrectionData, br) + for attr in attrs + PSY.get_transformer_winding(attr) == target || continue + return _compute_icd_factor(br, attr, winding_num) + end + return one(Float64) +end + """Ybus branch entries for an arc in the wye model of a `ThreeWindingTransformer`.""" -function ybus_branch_entries(tp::ThreeWindingTransformerWinding) +function ybus_branch_entries(tp::ThreeWindingTransformerWinding, correction::Float64) br = get_transformer(tp) winding_number = get_winding_number(tp) if winding_number == 1 - Y_t = 1 / (PSY.get_r_primary(br) + PSY.get_x_primary(br) * 1im) + Y_t = 1 / ((PSY.get_r_primary(br) + PSY.get_x_primary(br) * 1im) * correction) tap = PSY.get_primary_turns_ratio(br) * exp(PSY.get_α_primary(br) * 1im) c_tap = PSY.get_primary_turns_ratio(br) * exp(-1 * PSY.get_α_primary(br) * 1im) Y11 = Y_t / abs2(tap) @@ -472,9 +587,10 @@ function ybus_branch_entries(tp::ThreeWindingTransformerWinding) # primary bus alone gets the shunt term Y11 += y_shunt elseif winding_number == 2 - Y_t = 1 / (PSY.get_r_secondary(br) + PSY.get_x_secondary(br) * 1im) + Y_t = 1 / ((PSY.get_r_secondary(br) + PSY.get_x_secondary(br) * 1im) * correction) tap = PSY.get_secondary_turns_ratio(br) * exp(PSY.get_α_secondary(br) * 1im) - c_tap = PSY.get_secondary_turns_ratio(br) * exp(-1 * PSY.get_α_secondary(br) * 1im) + c_tap = + PSY.get_secondary_turns_ratio(br) * exp(-1 * PSY.get_α_secondary(br) * 1im) Y11 = Y_t / abs2(tap) if !isfinite(Y11) error( @@ -483,7 +599,7 @@ function ybus_branch_entries(tp::ThreeWindingTransformerWinding) ) end elseif winding_number == 3 - Y_t = 1 / (PSY.get_r_tertiary(br) + PSY.get_x_tertiary(br) * 1im) + Y_t = 1 / ((PSY.get_r_tertiary(br) + PSY.get_x_tertiary(br) * 1im) * correction) tap = PSY.get_tertiary_turns_ratio(br) * exp(PSY.get_α_tertiary(br) * 1im) c_tap = PSY.get_tertiary_turns_ratio(br) * exp(-1 * PSY.get_α_tertiary(br) * 1im) Y11 = Y_t / abs2(tap) @@ -500,6 +616,10 @@ function ybus_branch_entries(tp::ThreeWindingTransformerWinding) return (Y11, Y12, Y21, Y22) end +function ybus_branch_entries(tp::ThreeWindingTransformerWinding) + return ybus_branch_entries(tp, _get_impedance_correction_factor(tp)) +end + """Handles ybus entries for most 2-node AC branches. The types handled here are: `Line`, `DiscreteControlledACBranch`, `Transformer2W`, `TapTransformer`, and `PhaseShiftingTransformer`. """ @@ -509,6 +629,7 @@ function _ybus!( y21::Vector{YBUS_ELTYPE}, y22::Vector{YBUS_ELTYPE}, br::PSY.ACTransmission, + correction::Float64, num_bus::Dict{Int, Int}, branch_ix::Int, fb::Vector{Int}, @@ -521,12 +642,36 @@ function _ybus!( return end +function _ybus!( + y11::Vector{YBUS_ELTYPE}, + y12::Vector{YBUS_ELTYPE}, + y21::Vector{YBUS_ELTYPE}, + y22::Vector{YBUS_ELTYPE}, + br::PSY.TwoWindingTransformer, + correction::Float64, + num_bus::Dict{Int, Int}, + branch_ix::Int, + fb::Vector{Int}, + tb::Vector{Int}, + nr::NetworkReductionData, + adj::SparseArrays.SparseMatrixCSC{Int8, Int}, +) + add_branch_entries_to_indexing_maps!(num_bus, branch_ix, nr, fb, tb, adj, br) + Y11, Y12, Y21, Y22 = ybus_branch_entries(br, correction) + y11[branch_ix] = Y11 + y12[branch_ix] = Y12 + y21[branch_ix] = Y21 + y22[branch_ix] = Y22 + return +end + function _ybus!( y11::Vector{YBUS_ELTYPE}, y12::Vector{YBUS_ELTYPE}, y21::Vector{YBUS_ELTYPE}, y22::Vector{YBUS_ELTYPE}, br::PSY.DynamicBranch, + correction::Float64, num_bus::Dict{Int, Int}, branch_ix::Int, fb::Vector{Int}, @@ -534,19 +679,7 @@ function _ybus!( nr::NetworkReductionData, adj::SparseArrays.SparseMatrixCSC{Int8, Int}, ) - _ybus!( - y11, - y12, - y21, - y22, - br.branch, - num_bus, - branch_ix, - fb, - tb, - nr, - adj, - ) + _ybus!(y11, y12, y21, y22, br.branch, correction, num_bus, branch_ix, fb, tb, nr, adj) return end @@ -563,6 +696,9 @@ function _ybus!( ix::Int, nr::NetworkReductionData, adj::SparseArrays.SparseMatrixCSC{Int8, Int}, + c1::Float64, + c2::Float64, + c3::Float64, ) primary_star_arc = PSY.get_primary_star_arc(br) secondary_star_arc = PSY.get_secondary_star_arc(br) @@ -578,7 +714,8 @@ function _ybus!( adj[star_ix, primary_ix] = -1 fb[offset_ix + ix + n_entries] = primary_ix tb[offset_ix + ix + n_entries] = star_ix - (Y11, Y12, Y21, Y22) = ybus_branch_entries(ThreeWindingTransformerWinding(br, 1)) + (Y11, Y12, Y21, Y22) = + ybus_branch_entries(ThreeWindingTransformerWinding(br, 1), c1) y11[offset_ix + ix + n_entries] = Y11 y12[offset_ix + ix + n_entries] = Y12 y21[offset_ix + ix + n_entries] = Y21 @@ -591,7 +728,8 @@ function _ybus!( adj[star_ix, secondary_ix] = -1 fb[offset_ix + ix + n_entries] = secondary_ix tb[offset_ix + ix + n_entries] = star_ix - (Y11, Y12, Y21, Y22) = ybus_branch_entries(ThreeWindingTransformerWinding(br, 2)) + (Y11, Y12, Y21, Y22) = + ybus_branch_entries(ThreeWindingTransformerWinding(br, 2), c2) y11[offset_ix + ix + n_entries] = Y11 y12[offset_ix + ix + n_entries] = Y12 y21[offset_ix + ix + n_entries] = Y21 @@ -604,7 +742,8 @@ function _ybus!( adj[star_ix, tertiary_ix] = -1 fb[offset_ix + ix + n_entries] = tertiary_ix tb[offset_ix + ix + n_entries] = star_ix - (Y11, Y12, Y21, Y22) = ybus_branch_entries(ThreeWindingTransformerWinding(br, 3)) + (Y11, Y12, Y21, Y22) = + ybus_branch_entries(ThreeWindingTransformerWinding(br, 3), c3) y11[offset_ix + ix + n_entries] = Y11 y12[offset_ix + ix + n_entries] = Y12 y21[offset_ix + ix + n_entries] = Y21 @@ -678,6 +817,7 @@ function _buildybus!( fixed_admittances::Vector{PSY.FixedAdmittance}, switched_admittances::Vector{PSY.SwitchedAdmittance}, standard_loads::Vector{PSY.StandardLoad}, + icd_map::Union{Nothing, ICDCorrectionMap} = nothing, ) branch_entries_transformer_3w = 0 for br in transformer_3w @@ -701,11 +841,24 @@ function _buildybus!( y22 = zeros(YBUS_ELTYPE, branchcount) ysh = zeros(YBUS_ELTYPE, fa_count + sa_count + sl_count) - _foreach_ybus_branch(branches) do b, ix + _foreach_ybus_branch(branches, icd_map) do b, ix, correction if PSY.get_name(b) == "init" throw(DataFormatError("The data in Branch is invalid")) end - _ybus!(y11, y12, y21, y22, b, num_bus, ix, fb, tb, network_reduction_data, adj) + _ybus!( + y11, + y12, + y21, + y22, + b, + correction, + num_bus, + ix, + fb, + tb, + network_reduction_data, + adj, + ) end ix = 1 @@ -713,6 +866,14 @@ function _buildybus!( if PSY.get_name(b) == "init" throw(DataFormatError("The data in Transformer3W is invalid")) end + uuid = IS.get_uuid(b) + c1, c2, c3 = if isnothing(icd_map) + 1.0, 1.0, 1.0 + else + get(icd_map.map_3w, (uuid, 1), 1.0), + get(icd_map.map_3w, (uuid, 2), 1.0), + get(icd_map.map_3w, (uuid, 3), 1.0) + end n_entries = _ybus!( y11, y12, @@ -726,6 +887,9 @@ function _buildybus!( ix, network_reduction_data, adj, + c1, + c2, + c3, ) ix += n_entries end @@ -892,6 +1056,7 @@ function Ybus( bus_lookup[b] = ix end adj = SparseArrays.spdiagm(ones(Int8, busnumber)) + icd_map = _build_icd_correction_map(sys) branches = _get_ybus_two_terminal_ac_branches(sys) append!(branches.breaker_switches, breaker_switches) transformer_3W = @@ -915,6 +1080,7 @@ function Ybus( fixed_admittances, switched_admittances, standard_loads, + icd_map, ) ybus = SparseArrays.sparse( [fb; fb; tb; tb; sb], # row indices diff --git a/src/YbusACBranches.jl b/src/YbusACBranches.jl index 473ffd357..c5431c476 100644 --- a/src/YbusACBranches.jl +++ b/src/YbusACBranches.jl @@ -1,3 +1,8 @@ +struct ICDCorrectionMap + map_2w::Dict{Base.UUID, Float64} + map_3w::Dict{Tuple{Base.UUID, Int}, Float64} +end + struct YbusACBranches lines::Vector{PSY.Line} monitored_lines::Vector{PSY.MonitoredLine} @@ -32,6 +37,14 @@ function _populate_ybus_branch_vector!( return end +_get_correction(::Nothing, ::PSY.ACTransmission) = 1.0 +_get_correction(::Nothing, ::PSY.DynamicBranch) = 1.0 +_get_correction(::ICDCorrectionMap, ::PSY.ACTransmission) = 1.0 +_get_correction(map::ICDCorrectionMap, br::PSY.TwoWindingTransformer) = + get(map.map_2w, IS.get_uuid(br), 1.0) +_get_correction(map::ICDCorrectionMap, br::PSY.DynamicBranch) = + get(map.map_2w, IS.get_uuid(br.branch), 1.0) + function _get_ybus_two_terminal_ac_branches(sys::PSY.System)::YbusACBranches branches = YbusACBranches( Vector{PSY.Line}(), @@ -56,25 +69,27 @@ end function _foreach_ybus_branch( f::F, branches::YbusACBranches, + icd_map::Union{Nothing, ICDCorrectionMap}, ) where {F <: Function} - ix = _foreach_typed_branches(f, branches.lines, 0) - ix = _foreach_typed_branches(f, branches.monitored_lines, ix) - ix = _foreach_typed_branches(f, branches.generic_arc_impedances, ix) - ix = _foreach_typed_branches(f, branches.tap_transformers, ix) - ix = _foreach_typed_branches(f, branches.phase_shifting_transformers, ix) - ix = _foreach_typed_branches(f, branches.transformer2w, ix) - ix = _foreach_typed_branches(f, branches.dynamic_branches, ix) - ix = _foreach_typed_branches(f, branches.breaker_switches, ix) + ix = _foreach_typed_branches(f, branches.lines, icd_map, 0) + ix = _foreach_typed_branches(f, branches.monitored_lines, icd_map, ix) + ix = _foreach_typed_branches(f, branches.generic_arc_impedances, icd_map, ix) + ix = _foreach_typed_branches(f, branches.tap_transformers, icd_map, ix) + ix = _foreach_typed_branches(f, branches.phase_shifting_transformers, icd_map, ix) + ix = _foreach_typed_branches(f, branches.transformer2w, icd_map, ix) + ix = _foreach_typed_branches(f, branches.dynamic_branches, icd_map, ix) + ix = _foreach_typed_branches(f, branches.breaker_switches, icd_map, ix) return ix end function _foreach_typed_branches( f::F, vec::Vector{T}, + icd_map::Union{Nothing, ICDCorrectionMap}, offset::Int, ) where {F <: Function, T <: PSY.ACTransmission} for (i, br) in enumerate(vec) - f(br, offset + i) + f(br, offset + i, _get_correction(icd_map, br)) end return offset + length(vec) end diff --git a/test/test_data/modified_14bus_system_no_icd.txt b/test/test_data/modified_14bus_system_no_icd.txt new file mode 100644 index 000000000..84de40b1c --- /dev/null +++ b/test/test_data/modified_14bus_system_no_icd.txt @@ -0,0 +1,73 @@ + 107, 107, 3.10060620307922 , -10019.4659313190 + 107, 1100002, -0.00000000000000 , 10000.0000000000 + 107, 108, -3.10060620307922 , 19.4662113189697 + 1100002, 1100002, 0.00000000000000 , -30000.0000000000 + 1100002, 107, -0.00000000000000 , 10000.0000000000 + 1100002, 104, -0.00000000000000 , 10000.0000000000 + 1100002, 109, -0.00000000000000 , 10000.0000000000 + 103, 103, 132.904750823975 , -168.593001680914 + 103, 102, -19.0339088439941 , 120.448959350586 + 103, 104, -113.870841979980 , 48.1467323303223 + 102, 102, 5412.94931411743 , -4402.66704886226 + 102, 103, -19.0339088439941 , 120.448959350586 + 102, 104, -515.366577148438 , 376.782653808594 + 102, 101, -4878.04882812500 , 3902.43872070312 + 104, 104, 1028.59844207764 , -30440.8992626672 + 104, 1100002, -0.00000000000000 , 10000.0000000000 + 104, 103, -113.870841979980 , 48.1467323303223 + 104, 102, -515.366577148438 , 376.782653808594 + 104, 109, -0.00000000000000 , 10000.0000000000 + 104, 106, -0.00000000000000 , 10000.0000000000 + 104, 101, -399.361022949219 , 15.9744415283203 + 601, 601, 255.515943527222 , -266.450859436765 + 601, 301, -18.8712596893311 , 2.16857123374939 + 601, 109, -236.644683837891 , 264.337158203125 + 501, 501, 170.927175521851 , -68.8906958228908 + 501, 201, -157.686218261719 , 68.1742706298828 + 501, 109, -13.2409572601318 , 0.720955193042755 + 401, 401, 174.032964706421 , -239.601025610114 + 401, 200, -150.343994140625 , 102.884414672852 + 401, 109, -23.6889705657959 , 136.741210937500 + 301, 301, 24.4906101226807 , -51.2181056383997 + 301, 601, -18.8712596893311 , 2.16857123374939 + 301, 114, -5.61935043334961 , 49.1082344055176 + 201, 201, 253.106613159180 , -143.700978166540 + 201, 501, -157.686218261719 , 68.1742706298828 + 201, 114, -95.4203948974609 , 75.5306625366211 + 200, 200, 951.584594726562 , -167.455695092678 + 200, 401, -150.343994140625 , 102.884414672852 + 200, 114, -801.240600585938 , 64.6161804199219 + 114, 114, 902.280345916748 , -10189.1999023607 + 114, 301, -5.61935043334961 , 49.1082344055176 + 114, 201, -95.4203948974609 , 75.5306625366211 + 114, 200, -801.240600585938 , 64.6161804199219 + 114, 1100001, -0.00000000000000 , 10000.0000000000 + 112, 112, 76.9152069091797 , -10216.5844952832 + 112, 106, -76.9152069091797 , 216.588165283203 + 112, 1100001, -0.00000000000000 , 10000.0000000000 + 109, 109, 273.574611663818 , -30401.7677043328 + 109, 1100002, -0.00000000000000 , 10000.0000000000 + 109, 104, -0.00000000000000 , 10000.0000000000 + 109, 601, -236.644683837891 , 264.337158203125 + 109, 501, -13.2409572601318 , 0.720955193042755 + 109, 401, -23.6889705657959 , 136.741210937500 + 109, 110, -0.00000000000000 , 10000.0000000000 + 111, 111, 802.165511131287 , -149.270100883659 + 111, 110, -5.99848079681396 , 40.6076927185059 + 111, 106, -795.167030334473 , 110.667713165283 + 110, 110, 6.19848079979420 , -20036.6051827185 + 110, 109, -0.00000000000000 , 10000.0000000000 + 110, 111, -5.99848079681396 , 40.6076927185059 + 110, 1100001, -0.00000000000000 , 10000.0000000000 + 106, 106, 872.882237255573 , -10326.9494134365 + 106, 104, -0.00000000000000 , 10000.0000000000 + 106, 112, -76.9152069091797 , 216.588165283203 + 106, 111, -795.167030334473 , 110.667713165283 + 1100001, 1100001, 0.00000000000000 , -30000.0000000000 + 1100001, 114, -0.00000000000000 , 10000.0000000000 + 1100001, 112, -0.00000000000000 , 10000.0000000000 + 1100001, 110, -0.00000000000000 , 10000.0000000000 + + ZERO IMPEDANCE LINE CONNECTED BUSES (BUS IN MATRIX LISTED FIRST): + 104, 105 + 112, 113 diff --git a/test/test_data/modified_14bus_system_off_nominal_icd_ymatrix.txt b/test/test_data/modified_14bus_system_off_nominal_icd_ymatrix.txt new file mode 100644 index 000000000..603373836 --- /dev/null +++ b/test/test_data/modified_14bus_system_off_nominal_icd_ymatrix.txt @@ -0,0 +1,73 @@ + 107, 107, 3.10060620307922 , -10019.4659313190 + 107, 1100002, -0.00000000000000 , 10000.0000000000 + 107, 108, -3.10060620307922 , 19.4662113189697 + 1100002, 1100002, 0.00000000000000 , -30000.0000000000 + 1100002, 107, -0.00000000000000 , 10000.0000000000 + 1100002, 104, -0.00000000000000 , 10000.0000000000 + 1100002, 109, -0.00000000000000 , 10000.0000000000 + 103, 103, 132.904750823975 , -168.593001680914 + 103, 102, -19.0339088439941 , 120.448959350586 + 103, 104, -113.870841979980 , 48.1467323303223 + 102, 102, 5412.94931411743 , -4402.66704886226 + 102, 103, -19.0339088439941 , 120.448959350586 + 102, 104, -515.366577148438 , 376.782653808594 + 102, 101, -4878.04882812500 , 3902.43872070312 + 104, 104, 1028.59844207764 , -30149.6375439172 + 104, 1100002, -0.00000000000000 , 10000.0000000000 + 104, 103, -113.870841979980 , 48.1467323303223 + 104, 102, -515.366577148438 , 376.782653808594 + 104, 109, 0.00000000000000 , 10219.7246093750 + 104, 106, -0.00000000000000 , 10000.0000000000 + 104, 101, -399.361022949219 , 15.9744415283203 + 601, 601, 255.515943527222 , -266.450859436765 + 601, 301, -18.8712596893311 , 2.16857123374939 + 601, 109, -236.644683837891 , 264.337158203125 + 501, 501, 170.927175521851 , -68.8906958228908 + 501, 201, -157.686218261719 , 68.1742706298828 + 501, 109, -13.2409572601318 , 0.720955193042755 + 401, 401, 174.032964706421 , -239.601025610114 + 401, 200, -150.343994140625 , 102.884414672852 + 401, 109, -23.6889705657959 , 136.741210937500 + 301, 301, 24.4906101226807 , -51.2181056383997 + 301, 601, -18.8712596893311 , 2.16857123374939 + 301, 114, -5.61935043334961 , 49.1082344055176 + 201, 201, 253.106613159180 , -143.700978166540 + 201, 501, -157.686218261719 , 68.1742706298828 + 201, 114, -95.4203948974609 , 75.5306625366211 + 200, 200, 951.584594726562 , -167.455695092678 + 200, 401, -150.343994140625 , 102.884414672852 + 200, 114, -801.240600585938 , 64.6161804199219 + 114, 114, 902.280345916748 , -10189.1999023607 + 114, 301, -5.61935043334961 , 49.1082344055176 + 114, 201, -95.4203948974609 , 75.5306625366211 + 114, 200, -801.240600585938 , 64.6161804199219 + 114, 1100001, -0.00000000000000 , 10000.0000000000 + 112, 112, 76.9152069091797 , -9307.49367497070 + 112, 106, -76.9152069091797 , 216.588165283203 + 112, 1100001, -1578.61914062500 , 8952.79785156250 + 109, 109, 273.574611663818 , -30401.5167277703 + 109, 1100002, -0.00000000000000 , 10000.0000000000 + 109, 104, 0.00000000000000 , 10219.7246093750 + 109, 601, -236.644683837891 , 264.337158203125 + 109, 501, -13.2409572601318 , 0.720955193042755 + 109, 401, -23.6889705657959 , 136.741210937500 + 109, 110, 1953.09338378906 , 9033.41796875000 + 111, 111, 802.165511131287 , -149.270100883659 + 111, 110, -5.99848079681396 , 40.6076927185059 + 111, 106, -795.167030334473 , 110.667713165283 + 110, 110, 6.19848079979420 , -19278.7487374060 + 110, 109, -1953.09338378906 , 9033.41796875000 + 110, 111, -5.99848079681396 , 40.6076927185059 + 110, 1100001, -0.00000000000000 , 10000.0000000000 + 106, 106, 872.882237255573 , -10326.9494134365 + 106, 104, -0.00000000000000 , 10000.0000000000 + 106, 112, -76.9152069091797 , 216.588165283203 + 106, 111, -795.167030334473 , 110.667713165283 + 1100001, 1100001, 0.00000000000000 , -29090.9091796875 + 1100001, 114, -0.00000000000000 , 10000.0000000000 + 1100001, 112, 1578.61914062500 , 8952.79785156250 + 1100001, 110, -0.00000000000000 , 10000.0000000000 + + ZERO IMPEDANCE LINE CONNECTED BUSES (BUS IN MATRIX LISTED FIRST): + 104, 105 + 112, 113 diff --git a/test/test_ybus_psse.jl b/test/test_ybus_psse.jl index d9f1138a2..960ce1b53 100644 --- a/test/test_ybus_psse.jl +++ b/test/test_ybus_psse.jl @@ -85,7 +85,7 @@ end if col_bus ∈ keys(nr.reverse_bus_search_map) col_bus = nr.reverse_bus_search_map[col_bus] end - @test isapprox(Ybus_pnm[row_bus, col_bus], val, atol = 1e-3) + @test isapprox(Ybus_pnm[row_bus, col_bus], val, rtol = 2 * eps(Float32), atol = 0.0) end end @@ -118,7 +118,7 @@ end if col_bus ∈ keys(nr.reverse_bus_search_map) col_bus = nr.reverse_bus_search_map[col_bus] end - @test isapprox(Ybus_pnm[row_bus, col_bus], val, atol = 2e-3) + @test isapprox(Ybus_pnm[row_bus, col_bus], val, rtol = 4 * eps(Float32), atol = 0.0) end end @@ -202,7 +202,7 @@ end if col_bus ∈ keys(nr.reverse_bus_search_map) col_bus = nr.reverse_bus_search_map[col_bus] end - @test isapprox(Ybus_pnm[row_bus, col_bus], val; atol = 1e-1) + @test isapprox(Ybus_pnm[row_bus, col_bus], val; rtol = 9 * eps(Float32), atol = 0.0) end end @@ -238,6 +238,77 @@ end if col_bus ∈ keys(nr.reverse_bus_search_map) col_bus = nr.reverse_bus_search_map[col_bus] end - @test isapprox(Ybus_pnm[row_bus, col_bus], val; atol = 1e-5) + @test isapprox(Ybus_pnm[row_bus, col_bus], val; rtol = 2 * eps(Float32), atol = 0.0) + end +end + +@testset "psse_modified_14bus_off_nominal_icd" begin + # Three transformers have off-nominal settings that activate their ICTs: + # BUS 109-BUS 104-i_1 (TapTransformer): WINDV1=0.95, table 4 → F=1.030 + # BUS 110-BUS 109-i_1 (PhaseShiftingTransformer): ANG1=12.2°, table 3 → F=1.082 + # BUS 113-BUS 110-BUS 114-i_1 winding 1 (PST3W): ANG1=10.0°, table 9 → F=1.100 + sys = build_system(PSSEParsingTestSystems, "psse_modified_14bus_off_nominal_icd") + Ybus_pnm = Ybus(sys) + nr = Ybus_pnm.network_reduction_data + + Ybus_psse, b_ix_psse, row_buses, col_buses, y_values, reduced_bus_pairs_psse = + parse_psse_ybus( + joinpath(TEST_DATA_DIR, "modified_14bus_system_off_nominal_icd_ymatrix.txt"), + ) + for x in reduced_bus_pairs_psse + _test_psse_reduction_row(x, nr.reverse_bus_search_map) + end + + psse_ref_bus_numbers = [101, 108] + psse_skip_indices = indexin(psse_ref_bus_numbers, Ybus_pnm.axes[1]) + n_psse_excluded_elements = nnz(Ybus_pnm.data[psse_skip_indices, :]) + @test nnz(Ybus_pnm.data) == length(filter(!iszero, y_values)) + n_psse_excluded_elements + + psse_to_pnm_star_bus = Dict(1100001 => 1002, 1100002 => 1001) + for (row_bus, col_bus, val) in zip(row_buses, col_buses, y_values) + row_bus = get(psse_to_pnm_star_bus, row_bus, row_bus) + col_bus = get(psse_to_pnm_star_bus, col_bus, col_bus) + if row_bus ∈ keys(nr.reverse_bus_search_map) + row_bus = nr.reverse_bus_search_map[row_bus] + end + if col_bus ∈ keys(nr.reverse_bus_search_map) + col_bus = nr.reverse_bus_search_map[col_bus] + end + @test isapprox(Ybus_pnm[row_bus, col_bus], val; rtol = 2 * eps(Float32), atol = 0.0) + end +end + +@testset "psse_modified_14bus_nominal_icd" begin + sys = build_system(PSSEParsingTestSystems, "psse_modified_14bus_nominal_icd") + Ybus_pnm = Ybus(sys) + nr = Ybus_pnm.network_reduction_data + + Ybus_psse, b_ix_psse, row_buses, col_buses, y_values, reduced_bus_pairs_psse = + parse_psse_ybus(joinpath(TEST_DATA_DIR, "modified_14bus_system_no_icd.txt")) + for x in reduced_bus_pairs_psse + _test_psse_reduction_row(x, nr.reverse_bus_search_map) + end + + # PSS/E excluded rows for both type-3 buses it saw (101 and 108); count those + # PNM entries to reconcile the NNZ check. + psse_ref_bus_numbers = [101, 108] + psse_skip_indices = indexin(psse_ref_bus_numbers, Ybus_pnm.axes[1]) + n_psse_excluded_elements = nnz(Ybus_pnm.data[psse_skip_indices, :]) + @test nnz(Ybus_pnm.data) == length(filter(!iszero, y_values)) + n_psse_excluded_elements + + # PSS/E numbers its switch star-buses as 1100001/1100002; PNM assigns 1002/1001. + # Adjacency-derived mapping: PSS/E 1100001 ↔ {114,112,110} = PNM 1002; + # PSS/E 1100002 ↔ {107,104,109} = PNM 1001. + psse_to_pnm_star_bus = Dict(1100001 => 1002, 1100002 => 1001) + for (row_bus, col_bus, val) in zip(row_buses, col_buses, y_values) + row_bus = get(psse_to_pnm_star_bus, row_bus, row_bus) + col_bus = get(psse_to_pnm_star_bus, col_bus, col_bus) + if row_bus ∈ keys(nr.reverse_bus_search_map) + row_bus = nr.reverse_bus_search_map[row_bus] + end + if col_bus ∈ keys(nr.reverse_bus_search_map) + col_bus = nr.reverse_bus_search_map[col_bus] + end + @test isapprox(Ybus_pnm[row_bus, col_bus], val; rtol = 2 * eps(Float32), atol = 0.0) end end