ADD: center tap equivalence (#447)

* Add center tap equiv

* FIX: Julia 1.6 compat

* FIX: transformer test deepcopy

* DOC: changelog

---------

Co-authored-by: David M Fobes <dfobes@lanl.gov>

Author: jtabarez

CHANGELOG.md

@@ -2,7 +2,7 @@
 
 ## staged
 
-- none
+- Added transformation to remove center-tap (distribution) transformers and replace them with equivalenced lines
 
 ## v0.15.0
 

src/data_model/transformations/misc.jl

@@ -784,3 +784,102 @@ function join_lines!(data_eng::Dict{String,Any})
 
     return data_eng
 end
+
+
+"""
+    remove_distribution_transformers!(data_eng::Dict{String,<:Any})
+
+Removes distribution transformers and replaces them with lines. The current version only
+supports distribution transformers supplying loads i.e. no lines between transformer and load.
+"""
+function remove_distribution_transformers!(data_eng::Dict{String,<:Any})
+    if haskey(data_eng, "transformer")
+        for (_, transformer) in data_eng["transformer"]
+            buses = transformer["bus"]
+            line_connections = [0 for i = 1:length(buses)]
+            for (_, load) in data_eng["load"]
+                if load["bus"] in buses
+                    indxs = findall(x->x==load["bus"], buses)
+                    for indx in indxs
+                        indx != 0 ? line_connections[indx] = 1 : nothing
+                    end
+                end
+            end
+            for (_, line) in data_eng["line"]
+                if line["f_bus"] in buses
+                    indxs = findall(x->x==line["f_bus"], buses)
+                    for indx in indxs
+                        indx != 0 ? line_connections[indx] = 0 : nothing
+                    end
+                end
+                if line["t_bus"] in buses
+                    indxs = findall(x->x==line["t_bus"], buses)
+                    for indx in indxs
+                        indx != 0 ? line_connections[indx] = 0 : nothing
+                    end
+                end
+            end
+            if sum(line_connections) > 0
+                if length(buses) == 3
+                    _equivalance_center_tap!(transformer, data_eng)
+                elseif length(buses) == 2
+                    nothing
+                end
+            end
+        end
+    end
+end
+
+
+"""
+    _equivalance_center_tap!(transformer::Dict{String,<:Any}, data_eng::Dict{String,<:Any})
+
+Removes center tap transformers based on Eq. (1) from Kersting's paper
+'Center-Tapped Transformers and 120/240-V Secondary Models'
+Z0 = 0.5*r_t + j0.8*x_t
+"""
+function _equivalance_center_tap!(transformer, data_eng)
+    name = string("_virtual.", transformer["name"])
+    xs = transformer["xsc"][1] / 0.8
+    rw = transformer["rw"][1] / 0.5
+    data_eng["line"][name] = Dict{String,Any}(
+        "status" => deepcopy(transformer["status"]),
+        "length" => 1.0,
+        "units" => "ft",
+        "source_id" => transformer["source_id"],
+        "t_connections" => transformer["connections"][1],
+        "f_connections" => transformer["connections"][1],
+        "f_bus" => transformer["bus"][1],
+        "t_bus" => transformer["bus"][2],
+        "name" => name,
+        "rs" => fill(rw, 1, 1),
+        "xs" => fill(xs, 1, 1),
+        "g_fr" => fill(0.0, 1, 1),
+        "g_to" => fill(0.0, 1, 1),
+        "b_fr" => fill(0.0, 1, 1),
+        "b_to" => fill(0.0, 1, 1),
+    )
+    buses = Dict()
+    for bus in transformer["bus"]
+        buses[bus] = Dict("name" => string("_virtual.", bus), "replace" => false)
+    end
+    for (id, load) in data_eng["load"]
+        if load["bus"] in keys(buses)
+            buses[load["bus"]]["replace"] = true
+            _virtual_load = deepcopy(load)
+            _virtual_load["connections"] = transformer["connections"][1]
+            _virtual_load["vm_nom"] = transformer["vm_nom"][1]
+            _virtual_load["configuration"] = WYE
+            _virtual_load["pd_nom"] = [sum(_virtual_load["pd_nom"])]
+            _virtual_load["qd_nom"] = [sum(_virtual_load["qd_nom"])]
+            data_eng["load"][string("_virtual.", id)] = _virtual_load
+            load["status"] = DISABLED
+        end
+    end
+    for (id, bus) in buses
+        if bus["replace"]
+            data_eng["bus"][id]["terminals"] = transformer["connections"][1]
+        end
+    end
+    data_eng["transformer"][transformer["name"]]["status"] = DISABLED
+end

test/data/opendss/dist_transformer.dss

@@ -0,0 +1,67 @@
+clear
+
+! Base Frequency
+Set DefaultBaseFrequency=60
+
+! New Circuit
+New circuit.ut_trans
+~ BasekV=11 BaseMVA=0.1 pu=1.0  ISC3=9999999999 ISC1=9999999999 
+
+! Line Codes
+New LineCode.OH_3p nphases= 3 Units= mi
+ ~ Rmatrix= (1.191602 | 0.234849 1.191602 | 0.234849 0.234849 1.191602 ) 
+ ~ Xmatrix= (1.209452 | 0.489263 1.209452 | 0.489263 0.489263 1.209452 ) 
+New LineCode.OH_1p nphases= 1 Units= mi
+ ~ Rmatrix= (1.194260 ) 
+ ~ Xmatrix= (1.205420 ) 
+New LineCode.UG_3p nphases= 3 Units= mi
+ ~ Rmatrix= (1.009423 | 0.409732 1.009423 | 0.409732 0.409732 1.009423 ) 
+ ~ Xmatrix= (0.493164 | 0.100849 0.493164 | 0.100849 0.100849 0.493164 ) 
+ ~ Cmatrix= (286.101593 | 0.000000 286.101593 | 0.000000 0.000000 286.101593 )     
+
+! Lines
+New Line.LINE1 Bus1=SourceBus Bus2=1 phases=3 
+~ length=10 units=Mi LineCode=OH_3p
+New Line.LINE2 Bus1=1 Bus2=2 phases=3 
+~ length=5 units=Mi LineCode=UG_3p
+New Line.LINE3 Bus1=1 Bus2=3 phases=3 
+~ length=6 units=Mi LineCode=OH_3p
+New Line.LINE4 Bus1=3.1.0 Bus2=4.1.0 phases=1
+~ length=4 units=Mi LineCode=OH_1p
+New Line.LINE5 Bus1=3.2.0 Bus2=5.2.0 phases=1
+~ length=6 units=Mi Rmatrix= (1.194260 ) Xmatrix= (1.205420 ) 
+New Line.LINE6 Bus1=1 Bus2=6 phases=3
+~ length=5 units=Mi LineCode=OH_3p
+
+! Transformers 
+New Transformer.T1 Phases=1 Windings=3 XHL=2.256  XHT=2.256  XLT=1.504
+~ wdg=1 bus=4.1.0     conn=wye kV=6.35085    kva=25 %R=0.665
+~ wdg=2 bus=4_L.1.0 conn=wye kV=0.120     kva=25 %R=1.33 
+~ wdg=3 bus=4_L.0.2 conn=wye kV=0.120     kva=25 %R=1.33 
+New Transformer.T2 Phases=1 Windings=3 XHL=2.256  XHT=2.256  XLT=1.504
+~ wdg=1 bus=5.2.0     conn=wye kV=6.35085    kva=25 %R=0.665
+~ wdg=2 bus=5_L.1.0 conn=wye kV=0.120     kva=25 %R=1.33 
+~ wdg=3 bus=5_L.0.2 conn=wye kV=0.120     kva=25 %R=1.33
+New Transformer.T3 Phases=3 Windings=2 XHL=1.91
+~ wdg=1 bus=2.1.2.3.0     conn=wye kV=11  kva=75 %R=1.135
+~ wdg=2 bus=2_L.1.2.3.0 conn=wye kV=0.208 kva=75 %R=1.135
+New Transformer.T4 Phases=2 Windings=2 XHL=2.1
+~ wdg=1 bus=6.1.2     conn=delta kV=11  kva=75 %R=0.5
+~ wdg=2 bus=6_L.1.2 conn=delta kV=0.208 kva=75 %R=0.5
+
+! Load
+New  Load.L_6  phases=1  conn=delta  bus1=6_L.1.2  kV=0.208  kW=5.040000000000000    Kvar=2.296291839688011
+New  Load.L_2  phases=3  conn=wye  bus1=2_L.1.2.3.0  kV=0.208  kW=14.096000000000000    Kvar=6.422327335762342
+New  Load.L_4  phases=2  conn=wye  bus1=4_L.1.2  kV=0.208  kW=3.125000000000000    Kvar=1.027137828683947
+New  Load.L_5  phases=2  conn=wye  bus1=5_L.1.2  kV=0.208  kW=4.458000000000000    Kvar=1.195107182990382
+
+
+
+
+! Set Voltage Bases
+Set voltagebases=[11  .208]
+Calcvoltagebases
+
+! Solve network
+solve
+

test/test_cases.jl

@@ -61,3 +61,6 @@ test_gen_3ph_wye = parse_file("../test/data/en_validation_case_data/test_gen_3ph
 test_switch = parse_file("../test/data/en_validation_case_data/test_switch.dss", transformations=[remove_all_bounds!])
 test_gen_1ph_wye = parse_file("../test/data/en_validation_case_data/test_gen_1ph_wye.dss", transformations=[remove_all_bounds!])
 test_trans_dy = parse_file("../test/data/en_validation_case_data/test_trans_dy.dss", transformations=[remove_all_bounds!, transform_loops!])
+
+# distribution transformer equivalent cases
+dist_transformer = parse_file("../test/data/opendss/dist_transformer.dss")

test/transformer.jl

@@ -383,4 +383,29 @@
             @test all(isapprox.(result["solution"]["bus"]["tn_1"]["vm"], [1.045, 1.05]; atol=5E-3))
         end
     end
+
+    @testset "test center tap eq" begin
+        @testset "trans_3w_center_tap" begin
+            data = deepcopy(trans_3w_center_tap)
+            remove_distribution_transformers!(data)
+
+            @test data["transformer"]["xfmr_1"]["status"] == ENABLED
+            @test data["transformer"]["xfmr_2"]["status"] == ENABLED
+            @test data["transformer"]["xfmr_3"]["status"] == ENABLED
+        end
+
+        @testset "dist_transformer" begin
+            data = deepcopy(dist_transformer)
+            remove_distribution_transformers!(data)
+
+            result = solve_mc_pf(data, ACPUPowerModel, ipopt_solver; make_si=false)
+
+            @test data["transformer"]["t1"]["status"] == DISABLED
+            @test data["transformer"]["t2"]["status"] == DISABLED
+            @test norm(result["solution"]["bus"]["4"]["vm"]-[0.9990740842103211], Inf) <= 1.5E-5
+            @test norm(result["solution"]["bus"]["4"]["va"]-[-0.39064739635881085], Inf) <= 0.1
+            @test norm(result["solution"]["bus"]["4_l"]["vm"]-[0.9990723339621554], Inf) <= 1.5E-5
+            @test norm(result["solution"]["bus"]["4_l"]["va"]-[-0.3907533731198626], Inf) <= 0.1
+        end
+    end
 end