diff --git a/CHANGELOG.md b/CHANGELOG.md index c349ee822..b965560bc 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,11 @@ ## staged +- Added support for computing line constants from WireData, LineGeometry, LineSpacing, TSData and CNData +- Added Julia library SpecialFunctions for `besselj0` implementation +- Changed message that line is "being treated as superconducting" from `@info` to `@debug` +- Added support for WireData, LineGeometry, LineSpacing, TSData, and CNData dss objects +- Fixed bug in dss parser where when properties were assigned via `assign_property!`, the `prop_order` was not updated - Updated CI workflows to used Nodejs v16 scripts - Added UBF matrix power variables for switches [#423](https://github.com/lanl-ansi/PowerModelsDistribution.jl/issues/423) diff --git a/Project.toml b/Project.toml index 82f8dfd30..1499cf05e 100644 --- a/Project.toml +++ b/Project.toml @@ -16,18 +16,19 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Logging = "56ddb016-857b-54e1-b83d-db4d58db5568" LoggingExtras = "e6f89c97-d47a-5376-807f-9c37f3926c36" PolyhedralRelaxations = "2e741578-48fa-11ea-2d62-b52c946f73a0" +SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" [compat] CSV = "0.8.5, 0.9, 0.10" FilePaths = "0.8.3" Glob = "1.3" -InfrastructureModels = "0.7.3" -Ipopt = "0.9, 1.0.2" +InfrastructureModels = "0.7.3, 0.7.5" +Ipopt = "0.9, 1.0.2, 1.1" +LoggingExtras = "0.4.7, 1" JSON = "0.18, 0.19, 0.20, 0.21" JuMP = "0.22, 0.23, 1" -LoggingExtras = "0.4.7, 1" -PolyhedralRelaxations = "0.3.3" +PolyhedralRelaxations = "0.3.5" SCS = "0.9, 1.0, 1.1" julia = "1.6" diff --git a/src/PowerModelsDistribution.jl b/src/PowerModelsDistribution.jl index a3e4a972d..35fb1e1bb 100644 --- a/src/PowerModelsDistribution.jl +++ b/src/PowerModelsDistribution.jl @@ -13,6 +13,8 @@ module PowerModelsDistribution import InfrastructureModels + import SpecialFunctions + # Logging Utilities import Logging import LoggingExtras @@ -79,6 +81,7 @@ module PowerModelsDistribution include("core/relaxation_scheme.jl") include("io/utils.jl") + include("io/dss/line_constants.jl") include("io/dss/dss_parse.jl") include("io/dss/dss_data_structs.jl") include("io/dss/dss_node_structs.jl") diff --git a/src/core/constraint_template.jl b/src/core/constraint_template.jl index 583f23364..2e057b3d5 100644 --- a/src/core/constraint_template.jl +++ b/src/core/constraint_template.jl @@ -474,7 +474,7 @@ function constraint_mc_ohms_yt_from(pm::AbstractUnbalancedPowerModel, i::Int; nw t_idx = (i, t_bus, f_bus) if all(all(isapprox.(branch[k], 0.0)) for k in ["br_r", "br_x", "g_fr", "g_to", "b_fr", "b_to"]) - @info "branch $(branch["source_id"]) being treated as superconducting (effective zero impedance)" + @debug "branch $(branch["source_id"]) being treated as superconducting (effective zero impedance)" if !haskey(con(pm, nw), :branch_flow) con(pm, nw)[:branch_flow] = Dict{Int,Vector{Vector{<:JuMP.ConstraintRef}}}() end @@ -502,7 +502,7 @@ function constraint_mc_ohms_yt_to(pm::AbstractUnbalancedPowerModel, i::Int; nw:: t_idx = (i, t_bus, f_bus) if all(all(isapprox.(branch[k], 0.0)) for k in ["br_r", "br_x", "g_fr", "g_to", "b_fr", "b_to"]) - @info "branch $(branch["source_id"]) being treated as superconducting (effective zero impedance)" + @debug "branch $(branch["source_id"]) being treated as superconducting (effective zero impedance)" if !haskey(con(pm, nw), :branch_flow) con(pm, nw)[:branch_flow] = Dict{Int,Vector{Vector{<:JuMP.ConstraintRef}}}() end diff --git a/src/core/types.jl b/src/core/types.jl index 1657ecff0..5d83a79ba 100644 --- a/src/core/types.jl +++ b/src/core/types.jl @@ -84,11 +84,12 @@ An Enum to describe whether an object is enabled or disabled An Enum to describe the type of capcontrol, e.g., kvar, voltage etc. """ -@enum CapControlType CAP_CURRENT CAP_VOLTAGE CAP_REACTIVE_POWER CAP_DISABLED +@enum CapControlType CAP_CURRENT CAP_VOLTAGE CAP_REACTIVE_POWER CAP_DISABLED CAP_TIME @doc "Capacitor control based on current" CAP_CURRENT @doc "Capacitor control based on voltage" CAP_VOLTAGE @doc "Capacitor control based on total reactive power (directional)" CAP_REACTIVE_POWER @doc "Capacitor control disabled" CAP_DISABLED +@doc "Capacitor control based on time" CAP_TIME "Collection of the built-in Enums for PowerModelsDistribution" const PowerModelsDistributionEnums = Union{DataModel,LoadModel,ShuntModel,SwitchState,ControlMode,ConnConfig,Dispatchable,Status,CapControlType} diff --git a/src/data_model/eng2math.jl b/src/data_model/eng2math.jl index 33944b521..3fb94c5e4 100644 --- a/src/data_model/eng2math.jl +++ b/src/data_model/eng2math.jl @@ -677,6 +677,14 @@ function _map_eng2math_shunt!(data_math::Dict{String,<:Any}, data_eng::Dict{Stri "t_bus" => data_math["switch"][elem_id]["t_bus"] ) end + elseif dss_obj_type == "capacitor" + elem_id = first(filter(x->x.second["source_id"] == replace(math_obj["controls"]["element"], "capacitor"=>"shunt"), data_math["shunt"])).first + math_obj["controls"]["element"] = Dict{String,Any}( + "type" => "shunt", + "index" => data_math["shunt"][elem_id]["index"], + "f_bus" => data_math["shunt"][elem_id]["shunt_bus"], + "t_bus" => data_math["shunt"][elem_id]["shunt_bus"] + ) else elem_id = first(filter(x->x.second["source_id"] == math_obj["controls"]["element"], data_math["transformer"])).first math_obj["controls"]["element"] = Dict{String,Any}( diff --git a/src/data_model/units.jl b/src/data_model/units.jl index b790d0fb4..6de58cce1 100644 --- a/src/data_model/units.jl +++ b/src/data_model/units.jl @@ -484,6 +484,8 @@ function _rebase_pu_shunt!(shunt::Dict{String,<:Any}, vbase::Real, sbase::Real, shunt["controls"]["vmin"] = shunt["controls"]["vmin"]*shunt["controls"]["ptratio"]/(vbase*voltage_scale_factor) shunt["controls"]["vmax"] = shunt["controls"]["vmax"]*shunt["controls"]["ptratio"]/(vbase*voltage_scale_factor) end + elseif shunt["controls"]["type"] == CAP_TIME + # do nothing else for (idx,val) in enumerate(shunt["controls"]["type"]) if shunt["controls"]["voltoverride"][idx] diff --git a/src/io/dss/dss2eng.jl b/src/io/dss/dss2eng.jl index f2a17d221..03e389207 100644 --- a/src/io/dss/dss2eng.jl +++ b/src/io/dss/dss2eng.jl @@ -480,15 +480,33 @@ function _dss2eng_line!(data_eng::Dict{String,<:Any}, data_dss::Dict{String,<:An eng_obj["cm_ub"] = fill(defaults["emergamps"], ncond) end - if any(haskey(dss_obj, key) && _is_after_linecode(dss_obj["prop_order"], key) for key in ["r0", "r1", "rg", "rmatrix"]) || !haskey(dss_obj, "linecode") + from_geometry = haskey(dss_obj, "geometry") && _is_after_linecode(dss_obj["prop_order"], "geometry") + from_spacing = haskey(dss_obj, "wires") && _is_after_linecode(dss_obj["prop_order"], "wires") + from_cncables = haskey(dss_obj, "cndata") && _is_after_linecode(dss_obj["prop_order"], "cncables") + from_tscables = haskey(dss_obj, "tsdata") && _is_after_linecode(dss_obj["prop_order"], "tscables") + if from_geometry || from_spacing || from_cncables || from_tscables + z, y = calculate_line_constants(data_dss, defaults) + + rs, xs = real(z), imag(z) + g, b = real(y), imag(y) + + eng_obj["rs"] = rs + eng_obj["xs"] = xs + eng_obj["b_fr"] = b ./ 2.0 + eng_obj["b_to"] = b ./ 2.0 + eng_obj["g_fr"] = g ./ 2.0 + eng_obj["g_to"] = g ./ 2.0 + end + + if any(haskey(dss_obj, key) && (_is_after_linecode(dss_obj["prop_order"], key) && _is_after(dss_obj["prop_order"], key, "geometry")) for key in ["r0", "r1", "rg", "rmatrix"]) || (!haskey(dss_obj, "linecode") && !haskey(dss_obj, "geometry") && !haskey(dss_obj, "wires")) eng_obj["rs"] = reshape(defaults["rmatrix"], nphases, nphases) end - if any(haskey(dss_obj, key) && _is_after_linecode(dss_obj["prop_order"], key) for key in ["x0", "x1", "xg", "xmatrix"]) || !haskey(dss_obj, "linecode") + if any(haskey(dss_obj, key) && _is_after_linecode(dss_obj["prop_order"], key) && _is_after(dss_obj["prop_order"], key, "geometry") for key in ["x0", "x1", "xg", "xmatrix"]) || (!haskey(dss_obj, "linecode") && !haskey(dss_obj, "geometry") && !haskey(dss_obj, "wires")) eng_obj["xs"] = reshape(defaults["xmatrix"], nphases, nphases) end - if any(haskey(dss_obj, key) && _is_after_linecode(dss_obj["prop_order"], key) for key in ["b0", "b1", "c0", "c1", "cmatrix"]) || !haskey(dss_obj, "linecode") + if any(haskey(dss_obj, key) && _is_after_linecode(dss_obj["prop_order"], key) && _is_after(dss_obj["prop_order"], key, "geometry") for key in ["b0", "b1", "c0", "c1", "cmatrix"]) || (!haskey(dss_obj, "linecode") && !haskey(dss_obj, "geometry") && !haskey(dss_obj, "wires")) eng_obj["b_fr"] = reshape(defaults["cmatrix"], nphases, nphases) ./ 2.0 eng_obj["b_to"] = reshape(defaults["cmatrix"], nphases, nphases) ./ 2.0 eng_obj["g_fr"] = fill(0.0, nphases, nphases) @@ -857,7 +875,7 @@ function _dss2eng_regcontrol!(data_eng::Dict{String,<:Any}, data_dss::Dict{Strin defaults = _apply_ordered_properties(_create_regcontrol(id; _to_kwargs(dss_obj)...), dss_obj) nrw = get(data_dss["transformer"]["$(dss_obj["transformer"])"],"windings",2) - nphases = data_dss["transformer"]["$(dss_obj["transformer"])"]["phases"] + nphases = get(data_dss["transformer"]["$(dss_obj["transformer"])"], "phases", 3) eng_obj = Dict{String,Any}( "vreg" => [[w == defaults["winding"] && p == defaults["ptphase"] ? defaults["vreg"] : 0.0 for p in 1:nphases] for w in 1:nrw], @@ -911,6 +929,9 @@ function _dss2eng_capcontrol!(data_eng::Dict{String,<:Any}, data_dss::Dict{Strin eng_obj["onsetting"] = [p == defaults["ctphase"] ? defaults["onsetting"] : 0.0 for p in 1:nphases] eng_obj["offsetting"] = [p == defaults["ctphase"] ? defaults["offsetting"] : 0.0 for p in 1:nphases] eng_obj["ctratio"] = [p == defaults["ctphase"] ? defaults["ctratio"] : 0.0 for p in 1:nphases] + elseif type == CAP_TIME + eng_obj["type"] = type + eng_obj["terminal"] = defaults["terminal"] end if defaults["voltoverride"] @@ -941,6 +962,11 @@ function _dss2eng_capcontrol!(data_eng::Dict{String,<:Any}, data_dss::Dict{Strin eng_obj["offsetting"][defaults["ctphase"]] = defaults["offsetting"] eng_obj["ctratio"][defaults["ptphase"]] = defaults["ctratio"] end + if type == CAP_TIME + eng_obj["type"] = type + eng_obj["terminal"] = defaults["terminal"] + end + if defaults["voltoverride"] eng_obj["voltoverride"][defaults["ptphase"]] = defaults["voltoverride"] eng_obj["ptratio"][defaults["ptphase"]] = defaults["ptratio"] diff --git a/src/io/dss/dss_data_structs.jl b/src/io/dss/dss_data_structs.jl index 956df1f27..7a1dd099b 100644 --- a/src/io/dss/dss_data_structs.jl +++ b/src/io/dss/dss_data_structs.jl @@ -90,7 +90,6 @@ function _create_linecode(name::String=""; kwargs...)::Dict{String,Any} end - "Transformer codes contain all of the same properties as a transformer except bus, buses, bank, xfmrcode" function _create_xfmrcode(name::String=""; kwargs...) windings = isempty(name) ? 3 : get(kwargs, :windings, 2) @@ -208,7 +207,6 @@ function _create_xfmrcode(name::String=""; kwargs...) end - """ Creates a Dict{String,Any} containing all expected properties for a LoadShape element. See OpenDSS documentation for valid fields and ways to specify @@ -435,3 +433,281 @@ function _create_regcontrol(name::String=""; kwargs...)::Dict{String,Any} "like" => get(kwargs, :like, ""), ) end + + +""" +Creates a Dict{String,Any} containing all of the properties of LineGeometry. See +OpenDSS documentation for valid fields and ways to specify the different +properties. +""" +function _create_linegeometry(name::String=""; kwargs...)::Dict{String,Any} + nconds = get(kwargs, :nconds, 3) + conds = filter(x->startswith(string(x), "cond"), keys(kwargs)) + + temp = Dict{String,Any}( + "wires" => get(kwargs, :wires, fill("", nconds)), + "cncables" => get(kwargs, :cncables, fill("", nconds)), + "tscables" => get(kwargs, :tscables, fill("", nconds)), + "fconds" => collect(1:nconds), + "xs" => rand(nconds).*1e-100, + "hs" => rand(nconds).*1e-100, + "unitss" => fill("ft", nconds), + ) + + for cond in conds + suffix = kwargs[cond] == 1 ? "" : "_$(Int(kwargs[cond]))" + for key in [:wire, :x, :h, :units, :cncable, :tscable] + subkey = Symbol(string(key, suffix)) + if haskey(kwargs, subkey) + temp[string(key, "s")][Int(kwargs[cond])] = kwargs[subkey] + end + end + end + + for k in ["wires", "tscables", "cncables"] + if all(isempty(v) for v in temp[k]) + temp[k] = String[] + end + end + + Dict{String,Any}( + "name" => name, + "nconds" => get(kwargs, :nconds, 3), + "nphases" => get(kwargs, :nphases, 3), + "cond" => get(kwargs, :cond, 1), + "wire" => get(kwargs, :wire, ""), + "x" => get(kwargs, :x, rand()*1e-100) * _convert_to_meters[get(kwargs, :units, "ft")], + "h" => get(kwargs, :h, rand()*1e-100) * _convert_to_meters[get(kwargs, :units, "ft")], + "units" => "m", + "normamps" => get(kwargs, :normamps, 400.0), + "emergamps" => get(kwargs, :emergamps, 600.0), + "reduce" => get(kwargs, :reduce, false), + "spacing" => get(kwargs, :spacing, ""), + "wires" => temp["wires"], + "cncable" => get(kwargs, :cncable, ""), + "tscable" => get(kwargs, :tscable, ""), + "cncables" => temp["cncables"], + "tscables" => temp["tscables"], + "seasons" => get(kwargs, :seasons, 1), + "ratings" => get(kwargs, :ratings, Float64[600]), + "like" => get(kwargs, :like, ""), + # private properties + "fconds" => temp["fconds"], + "fx" => temp["xs"] .* [_convert_to_meters[u] for u in temp["unitss"]], + "fh" => temp["hs"] .* [_convert_to_meters[u] for u in temp["unitss"]], + "funits" => fill("m", nconds), + ) +end + + +""" +Creates a Dict{String,Any} containing all of the properties of WireData. See +OpenDSS documentation for valid fields and ways to specify the different +properties. +""" +function _create_wiredata(name::String=""; kwargs...)::Dict{String,Any} + if haskey(kwargs, :diam) || haskey(kwargs, :radius) + radunits = get(kwargs, :radunits, "none") + if haskey(kwargs, :diam) + diam = kwargs[:diam] + radius = kwargs[:diam] / 2 + elseif haskey(kwargs, :radius) + radius = kwargs[:radius] + diam = 2 * kwargs[:radius] + end + + if !haskey(kwargs, :gmrac) + gmrac = 0.7788 * kwargs[:radius] + gmrunits = get(kwargs, :radunits, "none") + else + gmrac = kwargs[:gmrac] + gmrunits = get(kwargs, :gmrunits, "none") + end + elseif haskey(kwargs, :gmrac) && !haskey(kwargs, :diag) && !haskey(kwargs, :radius) + radius = kwargs[:gmrac] / 0.7788 + diam = 2 * kwargs[:radius] + radunits = get(kwargs, :gmrunits, "none") + gmrac = kwargs[:gmrac] + gmrunits = radunits + else + radius = 1 + diam = 2 + gmrac = 0.7788 + radunits = "none" + gmrunits = "none" + end + + Dict{String,Any}( + "name" => name, + "rdc" => get(kwargs, :rdc, get(kwargs, :rac, 0.0)) / _convert_to_meters[get(kwargs, :runits, "none")], + "rac" => get(kwargs, :rac, get(kwargs, :rdc, 0.0)) / _convert_to_meters[get(kwargs, :runits, "none")], + "runits" => "m", + "gmrac" => gmrac * _convert_to_meters[gmrunits], + "gmrunits" => "m", + "radius" => radius * _convert_to_meters[radunits], + "capradius" => get(kwargs, :capradius, radius) * _convert_to_meters[radunits], + "radunits" => "m", + "normamps" => get(kwargs, :normamps, get(kwargs, :emergamps, 600.0) / 1.5), + "emergamps" => get(kwargs, :emergamps, get(kwargs, :normamps, 400.0) * 1.5), + "diam" => diam * _convert_to_meters[radunits], + "like" => get(kwargs, :like, ""), + ) +end + + +""" +Creates a Dict{String,Any} containing all of the properties of LineSpacing. See +OpenDSS documentation for valid fields and ways to specify the different +properties. +""" +function _create_linespacing(name::AbstractString=""; kwargs...)::Dict{String,Any} + Dict{String,Any}( + "name" => name, + "nconds" => get(kwargs, :nconds, 3), + "nphases" => get(kwargs, :nphases, get(kwargs, :nconds, 3)), + "x" => get(kwargs, :x, rand(get(kwargs, :nconds, 3)).*1e-100) .* _convert_to_meters[get(kwargs, :units, "ft")], + "h" => get(kwargs, :h, rand(get(kwargs, :nconds, 3)).*1e-100) .* _convert_to_meters[get(kwargs, :units, "ft")], + "units" => "m", + # internal properties + "fx" => get(kwargs, :x, rand(get(kwargs, :nconds, 3)).*1e-100) .* _convert_to_meters[get(kwargs, :units, "ft")], + "fh" => get(kwargs, :h, rand(get(kwargs, :nconds, 3)).*1e-100) .* _convert_to_meters[get(kwargs, :units, "ft")], + "funits" => fill("m", get(kwargs, :nconds, 3)) + ) +end + + +""" +Creates a Dict{String,Any} containing all of the properties of CNData. See +OpenDSS documentation for valid fields and ways to specify the different +properties. +""" +function _create_cndata(name::String=""; kwargs...)::Dict{String,Any} + if haskey(kwargs, :diam) || haskey(kwargs, :radius) + radunits = get(kwargs, :radunits, "none") + if haskey(kwargs, :diam) + diam = kwargs[:diam] + radius = kwargs[:diam] / 2 + elseif haskey(kwargs, :radius) + radius = kwargs[:radius] + diam = 2 * kwargs[:radius] + end + + if !haskey(kwargs, :gmrac) + gmrac = 0.7788 * kwargs[:radius] + gmrunits = get(kwargs, :radunits, "none") + else + gmrac = kwargs[:gmrac] + gmrunits = get(kwargs, :gmrunits, "none") + end + elseif haskey(kwargs, :gmrac) && !haskey(kwargs, :diam) && !haskey(kwargs, :radius) + radius = kwargs[:gmrac] / 0.7788 + diam = 2 * kwargs[:radius] + radunits = get(kwargs, :gmrunits, "none") + gmrac = kwargs[:gmrac] + gmrunits = radunits + else + radius = 1.0 + diam = 2.0 + gmrac = 0.7788 + radunits = "none" + gmrunits = "none" + end + + if haskey(kwargs, :diastrand) && !haskey(kwargs, :gmrstrand) + diastrand = kwargs[:diastrand] + gmrstrand = diastrand / 2 * 0.7788 + elseif haskey(kwargs, :gmrstrand) && !haskey(kwargs, :diastrand) + gmrstrand = kwargs[:gmrstrand] + diastrand = gmrstrand / 0.7788 * 2 + elseif haskey(kwargs, :diastrand) && haskey(kwargs, :gmrstrand) + diastrand = kwargs[:diastrand] + gmrstrand = kwargs[:gmrstrand] + else + diastrand = 2.0 + gmrstrand = 0.7788 + end + + Dict{String,Any}( + "name" => name, + "diacable" => get(kwargs, :diacable, 0.0) * _convert_to_meters[radunits], + "diains" => get(kwargs, :diains, 0.0) * _convert_to_meters[radunits], + "diam" => diam * _convert_to_meters[radunits], + "diastrand" => diastrand * _convert_to_meters[radunits], + "emergamps" => get(kwargs, :emergamps, get(kwargs, :normamps, 400.0) * 1.5), + "epsr" => get(kwargs, :epsr, 2.3), + "gmrac" => gmrac * _convert_to_meters[gmrunits], + "gmrstrand" => gmrstrand * _convert_to_meters[gmrunits], + "gmrunits" => "m", + "inslayer" => get(kwargs, :inslayer, 0.0) * _convert_to_meters[radunits], + "k" => get(kwargs, :k, 2), + "like" => get(kwargs, :like, ""), + "normamps" => get(kwargs, :normamps, get(kwargs, :emergamps, 600.0) / 1.5), + "rac" => get(kwargs, :rac, 1.02 * get(kwargs, :rdc, 1.0)) / _convert_to_meters[get(kwargs, :runits, "none")], + "radius" => radius * _convert_to_meters[radunits], + "radunits" => "m", + "rdc" => get(kwargs, :rdc, get(kwargs, :rac, 1.02) / 1.02) / _convert_to_meters[get(kwargs, :runits, "none")], + "rstrand" => get(kwargs, :rstrand, 0.0) / _convert_to_meters[get(kwargs, :runits, "none")], + "runits" => "m", + ) +end + + +""" +Creates a Dict{String,Any} containing all of the properties of TSData. See +OpenDSS documentation for valid fields and ways to specify the different +properties. +""" +function _create_tsdata(name::String=""; kwargs...)::Dict{String,Any} + if haskey(kwargs, :diam) || haskey(kwargs, :radius) + radunits = get(kwargs, :radunits, "none") + if haskey(kwargs, :diam) + diam = kwargs[:diam] + radius = kwargs[:diam] / 2 + elseif haskey(kwargs, :radius) + radius = kwargs[:radius] + diam = 2 * kwargs[:radius] + end + + if !haskey(kwargs, :gmrac) + gmrac = 0.7788 * kwargs[:radius] + gmrunits = get(kwargs, :radunits, "none") + else + gmrac = kwargs[:gmrac] + gmrunits = get(kwargs, :gmrunits, "none") + end + elseif haskey(kwargs, :gmrac) && !haskey(kwargs, :diag) && !haskey(kwargs, :radius) + radius = kwargs[:gmrac] / 0.7788 + diam = 2 * kwargs[:radius] + radunits = get(kwargs, :gmrunits, "none") + gmrac = kwargs[:gmrac] + gmrunits = radunits + else + radius = 1.0 + diam = 2.0 + gmrac = 0.7788 + radunits = "none" + gmrunits = "none" + end + + Dict{String,Any}( + "name" => name, + "diacable" => get(kwargs, :diacable, 0.0) * _convert_to_meters[radunits], + "diains" => get(kwargs, :diains, 0.0) * _convert_to_meters[radunits], + "diam" => diam * _convert_to_meters[radunits], + "diashield" => get(kwargs, :diashield, 0.0) * _convert_to_meters[radunits], + "emergamps" => get(kwargs, :emergamps, get(kwargs, :normamps, 400.0) * 1.5), + "epsr" => get(kwargs, :epsr, 2.3), + "gmrac" => gmrac * _convert_to_meters[gmrunits], + "gmrunits" => "m", + "inslayer" => get(kwargs, :inslayer, 0.0) * _convert_to_meters[radunits], + "like" => get(kwargs, :like, ""), + "normamps" => get(kwargs, :normamps, get(kwargs, :emergamps, 600.0) / 1.5), + "rac" => get(kwargs, :rac, get(kwargs, :rdc, 1.0) * 1.02) / _convert_to_meters[get(kwargs, :runits, "none")], + "radius" => radius * _convert_to_meters[radunits], + "radunits" => "m", + "rdc" => get(kwargs, :rdc, get(kwargs, :rac, 1.02) / 1.02) / _convert_to_meters[get(kwargs, :runits, "none")], + "runits" => "m", + "tapelap" => get(kwargs, :tapelap, 20.0), + "tapelayer" => get(kwargs, :tapelayer, 0.0) * _convert_to_meters[radunits], + ) +end diff --git a/src/io/dss/dss_parse.jl b/src/io/dss/dss_parse.jl index 90a3d6c23..f8cbeb2a0 100644 --- a/src/io/dss/dss_parse.jl +++ b/src/io/dss/dss_parse.jl @@ -278,6 +278,12 @@ const _dss_object_properties = Dict{String,Vector{String}}( ) +const _wdg_keys = Dict{String,Tuple{String,Vector{String}}}( + "transformer" => ("wdg", ["wdg", "bus", "conn", "kv", "kva", "tap", "%r", "rneut", "xneut"]), + "linegeometry" => ("cond", ["cond", "wire", "cncable", "tscable", "x", "h", "units"]) +) + + "parses single column load profile files" function _parse_csv_file(path::FilePaths.AbstractPath, type::AbstractString; header::Bool=false, column::Int=1, interval::Bool=false)::Union{Vector{String}, Tuple{Vector{String}, Vector{String}, Vector{String}}, Tuple{Vector{String}, Vector{String}}} open(first(Glob.glob([Glob.FilenameMatch(basename(path), "i")], dirname(path))), "r") do f @@ -637,21 +643,25 @@ the `key` and `value` of the property. If a property of the same name already exists inside `object`, the original value is converted to an array, and the new value is appended to the end. """ -function _add_property(object::Dict{String,<:Any}, key::SubString{String}, value::Any)::Dict{String,Any} +function _add_property(object::Dict{String,<:Any}, key::SubString{String}, value::Any; object_type::AbstractString="")::Dict{String,Any} if !haskey(object, "prop_order") object["prop_order"] = Vector{String}(["name"]) end - current_wdg = key == "wdg" ? value == "1" ? "" : "$value" : any(occursin("wdg", prop) for prop in object["prop_order"]) ? replace(split(filter(x->occursin("wdg", x), object["prop_order"])[end], "_")[end], "wdg"=>"") : "" - - if key in ["wdg", "bus", "conn", "kv", "kva", "tap", "%r", "rneut", "xneut"] - key = join(filter(p->!isempty(p), [key, current_wdg]), "_") + if haskey(_wdg_keys, object_type) + wdg_key, dup_keys = _wdg_keys[object_type] + current_wdg = key == wdg_key ? value == "1" ? "" : "$value" : any(occursin(wdg_key, prop) for prop in object["prop_order"]) ? replace(split(filter(x->occursin(wdg_key, x), object["prop_order"])[end], "_")[end], wdg_key=>"") : "" + if key in dup_keys + key = join(filter(p->!isempty(p), [key, current_wdg]), "_") + end end if haskey(object, lowercase(key)) - rmatch = match(r"_(\d+)$", key) - if typeof(rmatch) != Nothing - end_num = parse(Int, rmatch.captures[1]) + 1 + obj_keys = [k for k in keys(object) if startswith(k, key)] + rmatches = [match(r"_(\d+)$", k) for k in obj_keys] + ds = [parse(Int, rmatch.captures[1]) for rmatch in rmatches if typeof(rmatch) != Nothing] + if length(ds) > 0 + end_num = maximum(ds) + 1 key = replace(key, r"_(\d+)$" => "_$end_num") else key = string(key, "_2") @@ -722,7 +732,7 @@ function _parse_component(component::AbstractString, properties::AbstractString, value = _parse_mult_parameter(value; path=path) end - _add_property(object, key, value) + _add_property(object, key, value; object_type=obj_type) end return object @@ -788,6 +798,7 @@ of type `obj_type` named `obj_name` in `data_dss`. function _assign_property!(data_dss::Dict{String,<:Any}, obj_type::AbstractString, obj_name::AbstractString, property_name::AbstractString, property_value::Any) if haskey(data_dss, obj_type) && haskey(data_dss[obj_type], obj_name) data_dss[obj_type][obj_name][property_name] = property_value + push!(data_dss[obj_type][obj_name]["prop_order"], property_name) else @warn "Cannot find $obj_type object $obj_name." end diff --git a/src/io/dss/line_constants.jl b/src/io/dss/line_constants.jl new file mode 100644 index 000000000..39680cd8d --- /dev/null +++ b/src/io/dss/line_constants.jl @@ -0,0 +1,707 @@ +"Permeability of free space (N/A^2)" +const μ₀ = 4π*10^-7 + +"Permittivity of free space (F/m)" +const ε₀ = 8.8541878176e-12 + +"resistivity of copper tape shield (Ω-m)" +const ρₜₛ = 2.3718e-8 + + +"gets line geometry data for line, including applying line spacing if specified" +function _get_geometry_data(data_dss::Dict{String,<:Any}, geometry_id::String)::Dict{String,Any} + geom_obj = data_dss["linegeometry"][geometry_id] + _apply_like!(geom_obj, data_dss, "linegeometry") + + geometry = _apply_ordered_properties(_create_linegeometry(geometry_id; _to_kwargs(geom_obj)...), geom_obj) + + if !isempty(get(geometry, "spacing", "")) && !isempty(get(get(data_dss,"linespacing",Dict()), string(geometry["spacing"]), Dict())) + spacing = _get_spacing_data(data_dss, string(geometry["spacing"])) + + @assert geometry["nconds"] == spacing["nconds"] "Nconds on linegeometry.$(geometry_id) doesn't match nconds on linespacing.$(geometry["spacing"])" + + geometry["fx"] = spacing["fx"] + geometry["fh"] = spacing["fh"] + end + + return geometry +end + + +"get line spacing data for line or line geometry" +function _get_spacing_data(data_dss::Dict{String,<:Any}, spacing_id::String) + spacing_obj = data_dss["linespacing"][spacing_id] + _apply_like!(spacing_obj, data_dss, "linespacing") + + return _apply_ordered_properties(_create_linespacing(spacing_id; _to_kwargs(spacing_obj)...), spacing_obj) +end + + +"gets overhead wire data for line geometry" +function _get_wire_data(data_dss::Dict{String,<:Any}, wires::Vector{String})::Dict{String,Any} + wiredata = Dict{String,Any}(id => get(get(data_dss,"wiredata",Dict()), id, Dict{String,Any}()) for id in filter(x->!isempty(x), wires)) + @assert !isempty(wiredata) && all(!isempty(wd) for (_,wd) in wiredata) "Some wiredata is missing, cannot continue" + + return Dict{String,Any}( + id => _apply_ordered_properties(_create_wiredata(id; _to_kwargs(wd)...), wd) for (id,wd) in wiredata + ) +end + + +"gets concentric neutral cable data for line geometry" +function _get_cncable_data(data_dss::Dict{String,<:Any}, cncables::Vector{String})::Dict{String,Any} + cncabledata = Dict{String,Any}(id => get(get(data_dss,"cndata",Dict()), id, Dict{String,Any}()) for id in filter(x->!isempty(x), cncables)) + @assert !isempty(cncabledata) && all(!isempty(cncd) for (_,cncd) in cncabledata) "Some cndata is missing, cannot continue" + + return Dict{String,Any}( + id => _apply_ordered_properties(_create_cndata(id; _to_kwargs(cncd)...), cncd) for (id,cncd) in cncabledata + ) +end + + +"gets tape shielded cable data for line geometry" +function _get_tscable_data(data_dss::Dict{String,<:Any}, tscables::Vector{String})::Dict{String,Any} + tscabledata = Dict{String,Any}(id => get(get(data_dss,"tsdata",Dict()), id, Dict{String,Any}()) for id in filter(x->!isempty(x), tscables)) + @assert !isempty(tscabledata) && all(!isempty(tscd) for (_,tscd) in tscabledata) "Some tsdata is missing, cannot continue" + + return Dict{String,Any}( + id => _apply_ordered_properties(_create_tsdata(id; _to_kwargs(tscd)...), tscd) for (id,tscd) in tscabledata + ) +end + + +""" + calculate_line_constants(data_dss::Dict{String,<:Any}, line_defaults::Dict{String,<:Any})::Tuple{Matrix{Complex},Matrix{Complex}} + +Calculates line impedance and shunt admittance matrices for lines with line geometry, line spacing, wiredata, cncable, and/or tscable properties. +""" +function calculate_line_constants(data_dss::Dict{String,<:Any}, line_defaults::Dict{String,<:Any})::Tuple{Matrix{Complex},Matrix{Complex}} + geometry = !isempty(line_defaults["geometry"]) ? _get_geometry_data(data_dss, string(line_defaults["geometry"])) : missing + + cncables = !ismissing(geometry) && !isempty(geometry["cncables"]) ? filter(x->!isempty(x), geometry["cncables"]) : !isempty(line_defaults["cncables"]) ? filter(x->!isempty(x), line_defaults["cncables"]) : missing + ncncables = !ismissing(cncables) ? length(cncables) : missing + cndata = !ismissing(cncables) ? _get_cncable_data(data_dss, cncables) : missing + + tscables = !ismissing(geometry) && !isempty(geometry["tscables"]) ? filter(x->!isempty(x), geometry["tscables"]) : !isempty(line_defaults["tscables"]) ? filter(x->!isempty(x), line_defaults["tscables"]) : missing + ntscables = !ismissing(tscables) ? length(tscables) : missing + tsdata = !ismissing(tscables) ? _get_tscable_data(data_dss, tscables) : missing + + spacing = !isempty(line_defaults["spacing"]) ? _get_spacing_data(data_dss, string(line_defaults["spacing"])) : missing + + wires = !ismissing(geometry) && !isempty(geometry["wires"]) ? filter(x->!isempty(x), geometry["wires"]) : !isempty(line_defaults["wires"]) ? filter(x->!isempty(x), line_defaults["wires"]) : missing + nwires = !ismissing(wires) ? length(wires) : missing + wiredata = !ismissing(wires) ? _get_wire_data(data_dss, wires) : missing + + @assert count(.!(ismissing.([wires, cncables, tscables]))) > 0 "No wire, tscable, or cncable data is defined on line.$(line_defaults["name"])" + + ω = 2π * get(line_defaults, "basefreq", get(get(data_dss, "options", Dict()), "defaultbasefreq", 60.0)) + ω₀ = 2π * get(get(data_dss, "options", Dict()), "defaultbasefreq", 60.0) + + x = !ismissing(geometry) ? geometry["fx"] : !ismissing(spacing) ? spacing["fx"] : missing + y = !ismissing(geometry) ? geometry["fh"] : !ismissing(spacing) ? spacing["fh"] : missing + + gmr = !ismissing(wiredata) ? [wiredata[wire]["gmrac"] for wire in wires] : missing + capradius = !ismissing(wiredata) ? [wiredata[wire]["capradius"] for wire in wires] : missing + + rac = !ismissing(wiredata) ? [wiredata[wire]["rac"] for wire in wires] : missing + rdc = !ismissing(wiredata) ? [wiredata[wire]["rdc"] for wire in wires] : missing + + nphases = !ismissing(geometry) ? geometry["nphases"] : !ismissing(spacing) ? spacing["nphases"] : line_defaults["phases"] + nconds = !ismissing(geometry) ? geometry["nconds"] : !ismissing(spacing) ? spacing["nconds"] : nphases + + rac = !ismissing(cndata) ? [cndata[cable]["rac"] for cable in cncables] : rac + rdc = !ismissing(cndata) ? [cndata[cable]["rdc"] for cable in cncables] : rdc + gmr = !ismissing(cndata) ? [cndata[cable]["gmrac"] for cable in cncables] : gmr + + rstrand = !ismissing(cndata) ? [cndata[cable]["rstrand"] for cable in cncables] : missing + kstrand = !ismissing(cndata) ? [cndata[cable]["k"] for cable in cncables] : missing + diacable = !ismissing(cndata) ? [cndata[cable]["diacable"] for cable in cncables] : missing + diastrand = !ismissing(cndata) ? [cndata[cable]["diastrand"] for cable in cncables] : missing + gmrstrand = !ismissing(cndata) ? [cndata[cable]["gmrstrand"] for cable in cncables] : missing + epsr = !ismissing(cndata) ? [cndata[cable]["epsr"] for cable in cncables] : missing + diains = !ismissing(cndata) ? [cndata[cable]["diains"] for cable in cncables] : missing + inslayer = !ismissing(cndata) ? [cndata[cable]["inslayer"] for cable in cncables] : missing + + rac = !ismissing(tsdata) ? [tsdata[cable]["rac"] for cable in tscables] : rac + rdc = !ismissing(tsdata) ? [tsdata[cable]["rdc"] for cable in tscables] : rdc + gmr = !ismissing(tsdata) ? [tsdata[cable]["gmrac"] for cable in tscables] : gmr + + if count(.!(ismissing.([wires, cncables, tscables]))) > 1 + @assert sum(filter(x->!ismissing(x), [nwires, ncncables, ntscables])) == nconds "not enough wire/cable data for specified conductors on line.$(line_defaults["name"])" + push!(rac, [wiredata[wire]["rac"] for wire in wires]...) + push!(rdc, [wiredata[wire]["rdc"] for wire in wires]...) + push!(gmr, [wiredata[wire]["gmrac"] for wire in wires]...) + end + + diashield = !ismissing(tsdata) ? [tsdata[cable]["diashield"] for cable in tscables] : missing + tapelayer = !ismissing(tsdata) ? [tsdata[cable]["tapelayer"] for cable in tscables] : missing + tapelap = !ismissing(tsdata) ? [tsdata[cable]["tapelap"] for cable in tscables] : missing + + epsr = !ismissing(tsdata) ? [tsdata[cable]["epsr"] for cable in tscables] : epsr + diains = !ismissing(tsdata) ? [tsdata[cable]["diains"] for cable in tscables] : diains + inslayer = !ismissing(tsdata) ? [tsdata[cable]["inslayer"] for cable in tscables] : inslayer + + reduce = !ismissing(geometry) ? geometry["reduce"] : !ismissing(spacing) ? (line_defaults["phases"] == spacing["nconds"] ? false : true) : missing + + earth_model = !isempty(line_defaults["earthmodel"]) ? line_defaults["earthmodel"] : get(get(data_dss, "options", Dict()), "earthmodel", "deri") + + rho = line_defaults["rho"] + + Z, Y = calculate_line_constants( + x, + y, + ω, + gmr, + capradius, + nconds, + earth_model, + rac, + ω₀, + rdc, + rho, + nphases, + rstrand, + kstrand, + diacable, + diastrand, + gmrstrand, + epsr, + diains, + inslayer, + diashield, + tapelayer, + tapelap + ) + + if reduce + Z, Y = _kron(Z, Y, nphases) + end + + return Z, Y +end + + +""" + calculate_line_constants( + x::Vector{<:Real}, + y::Vector{<:Real}, + ω::Real, + gmr::Vector{<:Real}, + ::Union{Missing,Vector{<:Real}}, + nconds::Int, + earth_model::String, + R_ac::Vector{<:Real}, + ω₀::Real, + R_dc::Vector{<:Real}, + ρₑ::Real, + nphases::Int, + R_strand::Vector{<:Real}, + n_strand::Vector{<:Real}, + d_cable::Vector{<:Real}, + d_strand::Vector{<:Real}, + gmr_strand::Vector{<:Real}, + ε_ins::Vector{<:Real}, + d_ins::Vector{<:Real}, + t_ins::Vector{<:Real}, + ::Missing, + ::Missing, + ::Missing + )::Tuple{Matrix{Complex},Matrix{Complex}} + +Calculates the impedance and shunt admittance of concentric neutral cables + +# References + +- Nasser Tleis, Power Systems Modelling and Fault Analysis (Second Edition), Academic Press, 2019, ISBN 9780128151174. +- William H Kersting, Distribution System Modeling and Analysis (Forth Edition), CRC Press, 2018, ISBN 9781498772136. +- Andrea Ballanti, [Cable Modeling in OpenDSS](https://sourceforge.net/p/electricdss/code/HEAD/tree/trunk/Distrib/Doc/TechNote%20CableModelling.pdf?format=raw), 2017. +""" +function calculate_line_constants( + x::Vector{<:Real}, + y::Vector{<:Real}, + ω::Real, + gmr::Vector{<:Real}, + ::Union{Missing,Vector{<:Real}}, + nconds::Int, + earth_model::String, + R_ac::Vector{<:Real}, + ω₀::Real, + R_dc::Vector{<:Real}, + ρₑ::Real, + nphases::Int, + R_strand::Vector{<:Real}, + n_strand::Vector{<:Real}, + d_cable::Vector{<:Real}, + d_strand::Vector{<:Real}, + gmr_strand::Vector{<:Real}, + ε_ins::Vector{<:Real}, + d_ins::Vector{<:Real}, + t_ins::Vector{<:Real}, + ::Missing, + ::Missing, + ::Missing)::Tuple{Matrix{Complex},Matrix{Complex}} + + N = nconds + nphases + + Z = zeros(ComplexF64, N, N) + Y = zeros(ComplexF64, nphases, nphases) + + for i in 1:nconds + ii = i+nconds + for j in 1:nconds + jj = j+nconds + + Z_ije = Z_ije_ss = calc_earth_return_path_impedance(i, j, x, y, ρₑ, earth_model, ω, ω₀) + D_ij = sqrt((x[i]-x[j])^2+(y[i]-y[j])^2) + + if j <= i + if i == j + Z_ic = calc_internal_impedance(R_ac[i], R_dc[i], earth_model, ω).re + Z_ig = 1im * ω*μ₀/2π*log(1/gmr[i]) + + Z[i,i] = Z_ic + Z_ig + Z_ije + + if i <= nphases + r_cn = 1/2*(d_cable[i]-d_strand[i]) + gmr_cn = (gmr_strand[i]*n_strand[i]*(r_cn^(n_strand[i]-1)))^(1 / n_strand[i]) + + Z_ic_ss = R_strand[i] / n_strand[i] + Z_ig_ss = 1im * ω*μ₀/2π*log(1/gmr_cn) + + Z[ii,ii] = Z_ic_ss + Z_ig_ss + Z_ije_ss + end + else + Z_ijg = 1im * ω*μ₀/2π*log(1/D_ij) + + Z[i,j] = Z[j,i] = Z_ijg + Z_ije + + if i <= nphases + jj = j+nconds + + Z_ijg = 1im * ω*μ₀/2π*log(1/D_ij) + + Z[ii,jj] = Z[jj,ii] = Z_ijg + Z_ije + end + end + end + + if i <= nphases + r_cn = (d_cable[i]-d_strand[i])/2 + if i == j + D_ij = r_cn + else + D_ij = sqrt((x[i]-x[j])^2+(y[i]-y[j])^2) + D_ij = (D_ij^(n_strand[i]) - r_cn^(n_strand[i]))^(1/n_strand[i]) + end + + Z_ijg = 1im * ω*μ₀/2π*log(1/D_ij) + + Z[ii,j] = Z[j,ii] = Z_ijg + Z_ije + end + end + + if i <= nphases + r_outer = d_ins[i]/2 + r_inner = r_outer - t_ins[i] + Y[i,i] = 1im * 2π * ε₀ * ε_ins[i] * ω / log(r_outer / r_inner) + end + end + + Z = _kron(Z, nconds) + + C = Y .* 1e9 ./ ω + + return Z, C +end + + +""" + calculate_line_constants( + x::Vector{<:Real}, + y::Vector{<:Real}, + ω::Real, + gmr::Vector{<:Real}, + ::Union{Missing,Vector{<:Real}}, + nconds::Int, + earth_model::String, + R_ac::Vector{<:Real}, + ω₀::Real, + R_dc::Vector{<:Real}, + ρₑ::Real, + nphases::Int, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ε_ins::Vector{<:Real}, + d_ins::Vector{<:Real}, + t_ins::Vector{<:Real}, + d_shield::Vector{<:Real}, + t_tape::Vector{<:Real}, + lap_tape::Vector{<:Real} + )::Tuple{Matrix{Complex},Matrix{Complex}} + +Calculates the impedance and shunt admittance of tape shielded cables + +# References + +- Nasser Tleis, Power Systems Modelling and Fault Analysis (Second Edition), Academic Press, 2019, ISBN 9780128151174. +- William H Kersting, Distribution System Modeling and Analysis (Forth Edition), CRC Press, 2018, ISBN 9781498772136. +- Andrea Ballanti, [Cable Modeling in OpenDSS](https://sourceforge.net/p/electricdss/code/HEAD/tree/trunk/Distrib/Doc/TechNote%20CableModelling.pdf?format=raw), 2017. +""" +function calculate_line_constants( + x::Vector{<:Real}, + y::Vector{<:Real}, + ω::Real, + gmr::Vector{<:Real}, + ::Union{Missing,Vector{<:Real}}, + nconds::Int, + earth_model::String, + R_ac::Vector{<:Real}, + ω₀::Real, + R_dc::Vector{<:Real}, + ρₑ::Real, + nphases::Int, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ε_ins::Vector{<:Real}, + d_ins::Vector{<:Real}, + t_ins::Vector{<:Real}, + d_shield::Vector{<:Real}, + t_tape::Vector{<:Real}, + lap_tape::Vector{<:Real})::Tuple{Matrix{Complex},Matrix{Complex}} + + N = nconds + nphases + Z = zeros(Complex, N, N) + Y = zeros(Complex, nphases, nphases) + + for i in 1:nconds + ii = i+nconds + for j in 1:nconds + jj = j+nconds + + Z_ije = Z_ije_ss = calc_earth_return_path_impedance(i, j, x, y, ρₑ, earth_model, ω, ω₀) + D_ij = sqrt((x[i]-x[j])^2+(y[i]-y[j])^2) + + if j <= i + if i == j + Z_ic = calc_internal_impedance(R_ac[i], R_dc[i], earth_model, ω).re + Z_ig = 1im * ω*μ₀/2π*log(1/gmr[i]) + + Z[i,i] = Z_ic + Z_ig + Z_ije + + if i <= nphases + gmr_ts = (d_shield[i] - t_tape[i])/2 + + Z_ic_ss = (1+0im) * 0.3183 * ρₜₛ / (d_shield[i] * t_tape[i] * sqrt(50 / (100-lap_tape[i]))) + Z_ig_ss = (0+1im) * ω*μ₀/2π*log(1/gmr_ts) + + Z[ii,ii] = Z_ic_ss + Z_ig_ss + Z_ije_ss + end + else + Z_ijg = (0+1im) * ω*μ₀/2π*log(1/D_ij) + + Z[i,j] = Z[j,i] = Z_ijg + Z_ije + + if i <= nphases + jj = j+nconds + D_ij = sqrt((x[i]-x[j])^2+(y[i]-y[j])^2) + + Z_ijg = (0+1im) * ω*μ₀/2π*log(1 / D_ij) + + Z[ii,jj] = Z[jj,ii] = Z_ijg + Z_ije + end + end + end + + if i <= nphases + if i == j + D_ij = (d_shield[i] - t_tape[i])/2 + end + + Z_ijg = (0+1im) * ω*μ₀/2π*log(1/D_ij) + + Z[ii,j] = Z[j,ii] = Z_ijg + Z_ije + end + end + + if i <= nphases + r_outer = d_ins[i]/2 + r_inner = r_outer - t_ins[i] + Y[i,i] = (0+1im) * 2π * ε₀ * ε_ins[i] * ω / log(r_outer / r_inner) + end + end + + Z = _kron(Z, nconds) + + C = Y .* 1e9 ./ ω + + return Z, C +end + + +""" + calculate_line_constants( + x::Vector{<:Real}, + y::Vector{<:Real}, + ω::Real, + gmr::Vector{<:Real}, + r::Vector{<:Real}, + nconds::Int, + earth_model::String, + R_ac::Vector{<:Real}, + ω₀::Real, + R_dc::Vector{<:Real}, + ρₑ::Real, + ::Int, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing + )::Tuple{Matrix{Complex},Matrix{Complex}} + +Calculates impedance and shunt admittance for overhead lines. + +# References + +- Nasser Tleis, Power Systems Modelling and Fault Analysis (Second Edition), Academic Press, 2019, ISBN 9780128151174. +- William H Kersting, Distribution System Modeling and Analysis (Forth Edition), CRC Press, 2018, ISBN 9781498772136. +""" +function calculate_line_constants( + x::Vector{<:Real}, + y::Vector{<:Real}, + ω::Real, + gmr::Vector{<:Real}, + r::Vector{<:Real}, + nconds::Int, + earth_model::String, + R_ac::Vector{<:Real}, + ω₀::Real, + R_dc::Vector{<:Real}, + ρₑ::Real, + ::Int, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing, + ::Missing)::Tuple{Matrix{Complex},Matrix{Complex}} + + Z = zeros(ComplexF64, nconds, nconds) + P = zeros(ComplexF64, nconds, nconds) + + for i in 1:nconds + for j in 1:i + if i == j + Z_ic = calc_internal_impedance(R_ac[i], R_dc[i], earth_model, ω).re + # Tleis Eq 3.11 + Z_ig = 1im * ω*μ₀/2π*log(1/gmr[i]) + Z_ie = calc_earth_return_path_impedance(i,j, x, y, ρₑ, earth_model, ω, ω₀) + + # Tleis Eq 3.5b + Z[i,i] = Z_ic + Z_ig + Z_ie + + # Kersting Eq 5.9 + P[i,i] = 1im * -1/(2π * ε₀ * ω) * log(2 * y[i] / r[i]) + else + d_ij = sqrt((x[i]-x[j])^2+(y[i]-y[j])^2) + + # Tleis Eq 3.12b + X_ijg = (μ₀*ω/2π)*log(1/d_ij) + Z[i,j] = Z[j,i] = 1im * X_ijg + calc_earth_return_path_impedance(i,j, x, y, ρₑ, earth_model, ω, ω₀) + + # Kersting Pg 82 + S_ij = sqrt((x[i] - x[j])^2 + (y[i] + y[j])^2) + + # Kersting Eq 5.10 + P[i,j] = P[j,i] = 1im * -1/(2π * ε₀ * ω) * log(S_ij/d_ij) + end + end + end + + C = LinearAlgebra.pinv(P) .* 1e9 ./ ω + + return Z, C +end + + +""" + calc_earth_return_path_impedance(i::Int, j::Int, x::Vector{<:Real}, y::Vector{<:Real}, ρₑ::Real, earth_model::String, ω::Real, ω₀::Real) + +Calculates the earth return path impedance + +# References + +- Nasser Tleis, Power Systems Modelling and Fault Analysis (Second Edition), Academic Press, 2019, ISBN 9780128151174. +- A. Deri, G. Tevan, A. Semlyen, A. Castanheira, The Complex Ground Return Plane – A Simplified Model for Homogeneous and Multi-Layer Earth Return, IEEE Transactions on Power Apparatus and Systems, Vol. PAS-100, No. 8, pp. 3686-3693, August 1981. +""" +function calc_earth_return_path_impedance(i::Int, j::Int, x::Vector{<:Real}, y::Vector{<:Real}, ρₑ::Real, earth_model::String, ω::Real, ω₀::Real) + # Tleis Eq 3.7c + δ = 503.292*sqrt(ρₑ/(ω₀/2π)) + + # Tleis Eq 3.15 + D_erc = 1.309125*δ + + if earth_model ∈ ["simplecarson", "carson"] + Z_ij = (μ₀*ω/2π) * (π/4 + 1im * log(D_erc)) + elseif earth_model == "fullcarson" + # Tleis Pg 111 + b1 = 1/(3*sqrt(2)) + b2 = 1/16 + b3 = b1/(3*5) + b4 = b2/(4*6) + b5 = -b3/(5*7) + + c2 = 1.3659315 + c4 = c2 + 1/4 + 1/6 + + d2 = π/4*b2 + d4 = π/4*b4 + + + if i == j + # Tleis Eq 3.14a + D_ij = 2 * y[i] + # Tleis Eq 3.14b + θ_ij = 0 + else + # Tleis Eq 3.14a + D_ij = sqrt((x[i]-x[j])^2 + (y[i]+y[j])^2) + # Tleis Eq 3.14b + θ_ij = acos((y[i]+y[j])/D_ij) + end + + # Tleis Eq 3.14a + m_ij = sqrt(2) * D_ij / δ + + # Tleis Eq 3.13a + R_ije = (π/8 - b1*m_ij*cos(θ_ij) + b2*m_ij^2*(log(exp(c2)/m_ij)*cos(2*θ_ij) + θ_ij*sin(2*θ_ij)) + b3*m_ij^3*cos(3*θ_ij) - d4*m_ij^4*cos(4*θ_ij) - b5*m_ij^5*cos(5*θ_ij)) + + # Tleis Eq 3.13b + X_ije = (1/2*log(1.85138/m_ij) + b1*m_ij*cos(θ_ij) - d2*m_ij^2*cos(2*θ_ij) + b3*m_ij^3*cos(3*θ_ij) - b4*m_ij^4*(log(exp(c4)/m_ij)*cos(4*θ_ij)+θ_ij*sin(4*θ_ij)) + b5*m_ij^5*cos(5*θ_ij)) + + # Tleis Eq 3.17b + X_ije += 1/2*log(D_ij) + + Z_ij = ω/2π*μ₀ * 2 * (R_ije + 1im * X_ije) + elseif earth_model == "deri" + # Deri Eq 18 + p = 1 / sqrt(1im * ω * μ₀ / ρₑ) + + if i == j + # Deri Eq 3 + Z_ij = 1im * ω*μ₀/2π * log(2*(y[i] + p)) + else + # Deri Eq 4 + Z_ij = 1im * ω*μ₀/2π * log(sqrt((y[i] + y[j] + 2*p)^2 + (x[i]-x[j])^2)) + end + else + error("earth model $(earth_model) not recognized") + end + + return Z_ij +end + + +""" + calc_internal_impedance(R_ac::Real, R_dc::Real, earth_model::String, ω::Real) + +Calculates the internal impedance of the conductor + +# References + +- Nasser Tleis, Power Systems Modelling and Fault Analysis (Second Edition), Academic Press, 2019, ISBN 9780128151174. +""" +function calc_internal_impedance(R_ac::Real, R_dc::Real, earth_model::String, ω::Real) + if occursin("carson", earth_model) + L_i = μ₀/8π # internal inductance + Z_int = R_ac + 1im * ω*L_i + elseif earth_model == "deri" + δ = sqrt(R_dc*(ω/2π)*μ₀) + m = (1+1im)*sqrt(ω/2π*μ₀/R_dc) + Z_int = (1+1im) * (abs(m) > 35 ? (1+0im) : SpecialFunctions.besseli(0,m) / SpecialFunctions.besseli(1,m)) * δ / 2 + else + error("earth model $(earth_model) not recognized") + end + + return Z_int +end + + +""" + _kron(Z::Matrix{T}, Y::Matrix{T}, nconds::Int)::Tuple{Matrix{T}, Matrix{T}} where T <: Complex + +Kron reduces impedance and shunt admittance matrices down to size `(nconds, nconds)` +""" +function _kron(Z::Matrix{T}, Y::Matrix{T}, nconds::Int)::Tuple{Matrix{T}, Matrix{T}} where T <: Complex + Zout = deepcopy(Z) + Yout = deepcopy(Y) + + nrow, ncol = size(Z) + @assert nrow == ncol "cannot kron reduce non-square matrix" + + _Z = zeros(ComplexF64, nconds, nconds) + _Y = zeros(ComplexF64, nconds, nconds) + N = nrow + while N > nconds + _Z = zeros(ComplexF64, N-1, N-1) + _Y = zeros(ComplexF64, N-1, N-1) + + for i in 1:N-1 + for j in 1:N-1 + _Z[i,j] = Zout[i,j] - Zout[i,N]*Zout[j,N]/Zout[N,N] + _Y[i,j] = Yout[i,j] + end + end + + Zout = _Z + Yout = _Y + N -= 1 + end + + return Zout, Yout +end + + +""" + _kron(Z::Matrix{T}, nconds::Int)::Matrix{T} where T <: Complex + +Kron reduces impedance matrix down to size `(nconds, nconds)` +""" +function _kron(Z::Matrix{T}, nconds::Int)::Matrix{T} where T <: Complex + Zout = deepcopy(Z) + + nrow, ncol = size(Z) + @assert nrow == ncol "cannot kron reduce non-square matrix" + + _Z = zeros(T, nconds, nconds) + N = nrow + while N > nconds + _Z = zeros(T, N-1, N-1) + for i in 1:N-1 + for j in 1:N-1 + _Z[i,j] = Zout[i,j] - Zout[i,N]*Zout[j,N]/Zout[N,N] + end + end + Zout = _Z + N -= 1 + end + + return Zout +end diff --git a/src/io/utils.jl b/src/io/utils.jl index 1dab2aa89..27f9bb9b3 100644 --- a/src/io/utils.jl +++ b/src/io/utils.jl @@ -29,7 +29,8 @@ const _dss_supported_components = String[ "line", "linecode", "load", "generator", "capacitor", "reactor", "transformer", "pvsystem", "storage", "loadshape", "options", "xfmrcode", "vsource", "xycurve", "spectrum", "capcontrol", - "regcontrol", + "regcontrol", "linegeometry", "wiredata", "linespacing", + "cndata", "tsdata" ] "two number operators for reverse polish notation" @@ -91,6 +92,7 @@ const _dss2pmd_capcontrol_type = Dict{String,CapControlType}( "current" => CAP_CURRENT, "voltage" => CAP_VOLTAGE, ""=> CAP_DISABLED, + "time"=>CAP_TIME, ) diff --git a/test/data/opendss/IEEE13_Assets.dss b/test/data/opendss/IEEE13_Assets.dss new file mode 100644 index 000000000..c33dbf7e1 --- /dev/null +++ b/test/data/opendss/IEEE13_Assets.dss @@ -0,0 +1,142 @@ +// Reproduced from OpenDSS +// OpenDSS License +//* Copyright (c) 2008-2021, Electric Power Research Institute, Inc. +//* All rights reserved. +//* +//* Redistribution and use in source and binary forms, with or without +//* modification, are permitted provided that the following conditions are met: +//* * Redistributions of source code must retain the above copyright +//* notice, this list of conditions and the following disclaimer. +//* * Redistributions in binary form must reproduce the above copyright +//* notice, this list of conditions and the following disclaimer in the +//* documentation and/or other materials provided with the distribution. +//* * Neither the name of the Electric Power Research Institute, Inc., nor +//* the names of its contributors may be used to endorse or promote +//* products derived from this software without specific prior written +//* permission. +//* +//* THIS SOFTWARE IS PROVIDED BY Electric Power Research Institute, Inc., "AS IS" +//* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +//* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE +//* DISCLAIMED. IN NO EVENT SHALL Electric Power Research Institute, Inc., OR ANY OTHER +//* ENTITY CONTRIBUTING TO OR INVOLVED IN THE PROVISION OF THE SOFTWARE, BE +//* LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +//* CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +//* SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +//* INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +//* CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +//* ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +//* POSSIBILITY OF SUCH DAMAGE. + +Clear + +// added some CapControls and AssetInfosheets to test the CDPSM export + +new circuit.IEEE13NodecktAssets +~ basekv=115 pu=1.00 phases=3 bus1=SourceBus +~ Angle=30 ! advance angle 30 deg so result agree with published angle +~ MVAsc3=20000 MVASC1=21000 ! stiffen the source to approximate inf source + +New XfmrCode.SubXfmr Phases=3 Windings=2 XHL=(8 1000 /) +~ wdg=1 conn=delta kv=115 kva=5000 %r=(.5 1000 /) +~ wdg=2 conn=wye kv=4.16 kva=5000 %r=(.5 1000 /) +New XfmrCode.FdrXfmr Phases=3 Windings=2 XHL=2 %NoLoadLoss=0.6 %imag=1.1 +~ wdg=1 conn=Wye kv=4.16 kva=500 %r=.55 +~ wdg=2 conn=Wye kv=0.480 kva=500 %r=.55 +New XfmrCode.RegLeg phases=1 xhl=0.01 kvas=[1666 1666] kvs=[2.4 2.4] %LoadLoss=0.01 + +////////////////////////////////////////////////////////////////// + +new WireData.ACSR_556_5 NormAmps=730 DIAM=0.927 GMRac=0.37320 Rdc=0.035227273 Runits=kft Radunits=in gmrunits=in +new WireData.ACSR_4/0 NormAmps=340 DIAM=0.563 GMRac=0.09768 Rdc=0.112121212 Runits=kft Radunits=in gmrunits=in +new WireData.ACSR_1/0 NormAmps=230 DIAM=0.398 GMRac=0.05352 Rdc=0.212121212 Runits=kft Radunits=in gmrunits=in +new WireData.CU_1/0 NormAmps=100 DIAM=0.368 GMRac=0.13356 Rac=0.607 Runits=mi Radunits=in gmrunits=in +new TSData.TS_1/0 NormAmps=165 DIAM=0.368 GMRac=0.13320 Rac=0.97 Runits=mi Radunits=in gmrunits=in +~ EpsR=2.3 InsLayer=0.220 DiaIns=0.82 DiaCable=1.06 DiaShield=0.88 TapeLayer=0.005 TapeLap=20 +new CNData.CN_250 NormAmps=260 DIAM=0.567 GMRac=0.20520 Rac=0.41 Runits=mi Radunits=in gmrunits=in +~ EpsR=2.3 InsLayer=0.220 DiaIns=1.06 DiaCable=1.29 k=13 DiaStrand=0.0641 GmrStrand=0.02496 Rstrand=14.8722 + +new LineSpacing.500 nconds=4 nphases=3 units=ft x=[-4 -1 3 0] h=[28 28 28 24] +new LineSpacing.505 nconds=3 nphases=2 units=ft x=[-4 3 0] h=[28 28 24] +new LineSpacing.510 nconds=2 nphases=1 units=ft x=[0.5 0] h=[29 24] + +new LineGeometry.601 nconds=4 nphases=3 reduce=y spacing=500 wires=[ACSR_556_5 ACSR_556_5 ACSR_556_5 ACSR_4/0] +new LineGeometry.602 nconds=4 nphases=3 reduce=y spacing=500 wires=[ACSR_4/0 ACSR_4/0 ACSR_4/0 ACSR_4/0] +new LineGeometry.603 nconds=3 nphases=2 reduce=y spacing=505 wires=[ACSR_1/0 ACSR_1/0 ACSR_1/0] +new LineGeometry.604 like=603 +new LineGeometry.605 nconds=2 nphases=1 reduce=y spacing=510 wires=[ACSR_1/0 ACSR_1/0] + +new LineGeometry.606 nconds=3 nphases=3 reduce=y +~ cond=1 cncable=CN_250 x=-0.5 h=-4 units=ft +~ cond=2 cncable=CN_250 x= 0.0 h=-4 units=ft +~ cond=3 cncable=CN_250 x= 0.5 h=-4 units=ft + +new LineGeometry.607 nconds=2 nphases=1 reduce=y +~ cond=1 tscable=TS_1/0 x= 0.0000 h=-4 units=ft +~ cond=2 wire =CU_1/0 x= 0.2500 h=-4 units=ft + +////////////////////////////////////////////////////////////////// + +New Transformer.Sub XfmrCode=SubXfmr Buses=[SourceBus 650] +New Transformer.XFM1 XfmrCode=FdrXfmr Buses=[633 634] + +New Transformer.Reg1 XfmrCode=RegLeg Bank=Reg Buses=[650.1 RG60.1] +New Transformer.Reg2 XfmrCode=RegLeg Bank=Reg Buses=[650.2 RG60.2] +New Transformer.Reg3 XfmrCode=RegLeg Bank=Reg Buses=[650.3 RG60.3] + +new regcontrol.Reg1 transformer=Reg1 winding=2 vreg=122 band=2 ptratio=20 ctprim=700 R=3 X=9 +new regcontrol.Reg2 transformer=Reg2 winding=2 vreg=122 band=2 ptratio=20 ctprim=700 R=3 X=9 +new regcontrol.Reg3 transformer=Reg3 winding=2 vreg=122 band=2 ptratio=20 ctprim=700 R=3 X=9 + +New Load.671 Bus1=671.1.2.3 Phases=3 Conn=Delta Model=1 kV=4.16 kW=1155 kvar=660 +New Load.634a Bus1=634.1 Phases=1 Conn=Wye Model=1 kV=0.277 kW=160 kvar=110 +New Load.634b Bus1=634.2 Phases=1 Conn=Wye Model=1 kV=0.277 kW=120 kvar=90 +New Load.634c Bus1=634.3 Phases=1 Conn=Wye Model=1 kV=0.277 kW=120 kvar=90 +New Load.645 Bus1=645.2 Phases=1 Conn=Wye Model=1 kV=2.4 kW=170 kvar=125 +New Load.646 Bus1=646.2.3 Phases=1 Conn=Delta Model=2 kV=4.16 kW=230 kvar=132 +New Load.692 Bus1=692.3.1 Phases=1 Conn=Delta Model=5 kV=4.16 kW=170 kvar=151 +New Load.675a Bus1=675.1 Phases=1 Conn=Wye Model=1 kV=2.4 kW=485 kvar=190 +New Load.675b Bus1=675.2 Phases=1 Conn=Wye Model=1 kV=2.4 kW=68 kvar=60 +New Load.675c Bus1=675.3 Phases=1 Conn=Wye Model=1 kV=2.4 kW=290 kvar=212 +New Load.611 Bus1=611.3 Phases=1 Conn=Wye Model=5 kV=2.4 kW=170 kvar=80 +New Load.652 Bus1=652.1 Phases=1 Conn=Wye Model=2 kV=2.4 kW=128 kvar=86 +New Load.670a Bus1=670.1 Phases=1 Conn=Wye Model=1 kV=2.4 kW=17 kvar=10 +New Load.670b Bus1=670.2 Phases=1 Conn=Wye Model=1 kV=2.4 kW=66 kvar=38 +New Load.670c Bus1=670.3 Phases=1 Conn=Wye Model=1 kV=2.4 kW=117 kvar=68 + +!Bus 670 is the concentrated point load of the distributed load on line 632 to 671 located at 1/3 the distance from node 632 + +New Line.650632 Phases=3 Bus1=RG60.1.2.3 Bus2=632.1.2.3 Geometry=601 Length=2000 units=ft +New Line.632670 Phases=3 Bus1=632.1.2.3 Bus2=670.1.2.3 Geometry=601 Length=667 units=ft +New Line.670671 Phases=3 Bus1=670.1.2.3 Bus2=671.1.2.3 Geometry=601 Length=1333 units=ft +New Line.671680 Phases=3 Bus1=671.1.2.3 Bus2=680.1.2.3 Geometry=601 Length=1000 units=ft +New Line.632633 Phases=3 Bus1=632.1.2.3 Bus2=633.1.2.3 Geometry=602 Length=500 units=ft +New Line.632645 Phases=2 Bus1=632.3.2 Bus2=645.3.2 Geometry=603 Length=500 units=ft +// New Line.645646 Phases=2 Bus1=645.3.2 Bus2=646.3.2 Geometry=603 Length=300 units=ft +New Line.645646 Phases=2 Bus1=645.3.2 Bus2=646.3.2 Length=0.0568 units=mi normamps=230 +~ rmatrix = (1.3238 | 0.2066 1.3294 ) +~ xmatrix = (1.3569 | 0.4591 1.3471 ) +New Line.692675 Phases=3 Bus1=692.1.2.3 Bus2=675.1.2.3 Geometry=606 Length=500 units=ft +New Line.671684 Phases=2 Bus1=671.1.3 Bus2=684.1.3 Geometry=604 Length=300 units=ft +// New Line.684611 Phases=1 Bus1=684.3 Bus2=611.3 Geometry=605 Length=300 units=ft +New Line.684611 Phases=1 Bus1=684.3 Bus2=611.3 Spacing=510 Wires=[ACSR_1/0 ACSR_1/0] Length=300 units=ft +New Line.684652 Phases=1 Bus1=684.1 Bus2=652.1 Geometry=607 Length=800 units=ft +New Line.671692 Phases=3 Bus1=671 Bus2=692 Switch=y r1=1e-4 r0=1e-4 x1=0.000 x0=0.000 c1=0.000 c0=0.000 + +New Capacitor.Cap1 Bus1=675 phases=3 kVAR=600 kV=4.16 +New Capacitor.Cap2 Bus1=611.3 phases=1 kVAR=100 kV=2.4 +New CapControl.Cap1 capacitor=cap1 type=time on=8 off=19 element=capacitor.cap1 +New CapControl.Cap2 capacitor=cap2 type=voltage on=115 off=125 ptratio=20 ptphase=1 element=line.684611 terminal=2 + + +Set Voltagebases=[115, 4.16, .48] +calcv +Solve +//BusCoords IEEE13Node_BusXY.csv + +Transformer.Reg1.Taps=[1.0 1.0625] +Transformer.Reg2.Taps=[1.0 1.0500] +Transformer.Reg3.Taps=[1.0 1.06875] +Set Controlmode=OFF + +Solve diff --git a/test/data/opendss/test2_master.dss b/test/data/opendss/test2_master.dss index 46329c251..506d7faaf 100644 --- a/test/data/opendss/test2_master.dss +++ b/test/data/opendss/test2_master.dss @@ -87,7 +87,7 @@ new linespacing.test_spacing new wiredata.wire1 new wiredata.wire2 -new line.l9 bus1=b10 bus2=b11 spacing=test_spacing wires=['wire1', 'wire2'] length=10 +new line.l9 bus1=b10 bus2=b11 spacing=test_spacing wires=['wire1', 'wire2', 'wire1', 'wire2'] length=10 new storage.s1 bus1=b1 new pvsystem.solar1 bus1=b1 diff --git a/test/line_constants.jl b/test/line_constants.jl new file mode 100644 index 000000000..3439b49bd --- /dev/null +++ b/test/line_constants.jl @@ -0,0 +1,61 @@ +@testset "line constants" begin + eng = parse_file("../test/data/opendss/IEEE13_Assets.dss") + + z = [ + 0.000211485 9.80031e-5 9.61061e-5; + 9.80031e-5 0.000217924 9.93019e-5; + 9.61061e-5 9.93019e-5 0.000213946 + ] + 1im * [ + 0.000653882 0.000299058 0.000241464; + 0.000299058 0.000632756 0.000273149; + 0.000241464 0.000273149 0.000645754 + ] + y = [ + 0 0 0; + 0 0 0; + 0 0 0 + ] + 1im * [ + 0.00961082 -0.00291343 -0.00123058; + -0.00291343 0.0103002 -0.00230096; + -0.00123058 -0.00230096 0.00939874 + ] + + @test all(isapprox.(eng["line"]["650632"]["rs"] + 1im * eng["line"]["650632"]["xs"], z; atol=1e-4)) + @test all(isapprox.(eng["line"]["650632"]["g_fr"].*2 + 1im * eng["line"]["650632"]["b_fr"].*2, y; atol=1e-4)) + + z = [ + 0.000491885 0.000198743 0.000177491; + 0.000198743 0.000486187 0.000198743; + 0.000177491 0.000198743 0.000491885 + ] + 1im * [ + 0.00027695 1.99955e-5 -9.24851e-6; + 1.99955e-5 0.000250698 1.99955e-5; + -9.24851e-6 1.99955e-5 0.00027695 + ] + y = [ + 0 0 0; + 0 0 0; + 0 0 0 + ] + 1im * [ + 0.238581 0.0 0.0; + 0.0 0.238581 0.0; + 0.0 0.0 0.238581 + ] + + @test all(isapprox.(eng["line"]["692675"]["rs"] + 1im * eng["line"]["692675"]["xs"], z; atol=1e-5)) + @test all(isapprox.(eng["line"]["692675"]["g_fr"].*2 + 1im * eng["line"]["692675"]["b_fr"].*2, y; atol=1e-5)) + + z = [0.000802907;;] + 1im * [0.000429279;;] + y = [0;;] + 1im * [0.166359;;] + + @test all(isapprox.(eng["line"]["684652"]["rs"] + 1im * eng["line"]["684652"]["xs"], z; atol=1e-5)) + @test all(isapprox.(eng["line"]["684652"]["g_fr"].*2 + 1im * eng["line"]["684652"]["b_fr"].*2, y; atol=1e-5)) + + r = solve_mc_pf(eng, ACRUPowerModel, ipopt_solver; make_si=false, solution_processors=[sol_data_model!]) + + @test all(isapprox.(r["solution"]["bus"]["652"]["vm"], [0.974806]; atol=1e-4)) + @test all(isapprox.(r["solution"]["bus"]["652"]["va"], [-5.35]; atol=1e-0)) + + @test all(isapprox.(r["solution"]["bus"]["675"]["vm"], [0.975987, 1.04705, 0.984996]; atol=1e-4)) + @test all(isapprox.(r["solution"]["bus"]["675"]["va"], [-5.52, -121.71, 116.03]; atol=1e-0)) +end diff --git a/test/opendss.jl b/test/opendss.jl index bc5d6c8f5..06ed268c9 100644 --- a/test/opendss.jl +++ b/test/opendss.jl @@ -160,7 +160,7 @@ @testset "opendss parse line parsing wires - spacing properties" begin dss_data = parse_dss("../test/data/opendss/test2_master.dss") - @test isa(dss_data["line"]["l9"]["wires"], Vector{String}) && all(dss_data["line"]["l9"]["wires"] .== ["wire1", "wire2"]) + @test isa(dss_data["line"]["l9"]["wires"], Vector{String}) && all(dss_data["line"]["l9"]["wires"] .== ["wire1", "wire2", "wire1", "wire2"]) @test dss_data["line"]["l9"]["spacing"] == "test_spacing" @test haskey(dss_data, "wiredata") && (haskey(dss_data["wiredata"], "wire1") && haskey(dss_data["wiredata"], "wire2")) @test haskey(dss_data, "linespacing") && haskey(dss_data["linespacing"], "test_spacing") diff --git a/test/runtests.jl b/test/runtests.jl index 49950fc40..296fc8933 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -62,4 +62,6 @@ include("test_cases.jl") include("en_opf_bounds.jl") include("en_pf_validation.jl") + + include("line_constants.jl") end