e2nIEE / pandapower
1
# -*- coding: utf-8 -*-
2

3
# Copyright (c) 2016-2021 by University of Kassel and Fraunhofer Institute for Energy Economics
4
# and Energy System Technology (IEE), Kassel. All rights reserved.
5

6 1
import math
7 1
import numpy as np
8 1
import pandapower.auxiliary as aux
9 1
from pandapower.build_bus import _build_bus_ppc
10 1
from pandapower.build_gen import _build_gen_ppc
11
#from pandapower.pd2ppc import _ppc2ppci, _init_ppc
12 1
from pandapower.pypower.idx_brch import BR_B, BR_R, BR_X, F_BUS, T_BUS, branch_cols, BR_STATUS, SHIFT, TAP
13 1
from pandapower.pypower.idx_bus import BASE_KV, BS, GS
14 1
from pandapower.build_branch import _calc_tap_from_dataframe, _transformer_correction_factor, _calc_nominal_ratio_from_dataframe
15 1
from pandapower.build_branch import _switch_branches, _branches_with_oos_buses, _initialize_branch_lookup, _end_temperature_correction_factor
16

17 1
def _pd2ppc_zero(net, sequence=0):
18 1
    from pandapower.pd2ppc import _ppc2ppci, _init_ppc
19

20
    """
21
    Builds the ppc data structure for zero impedance system. Includes the impedance values of
22
    lines and transformers, but no load or generation data.
23

24
    For short-circuit calculation, the short-circuit impedance of external grids is also considered.
25
    """
26
    # select elements in service (time consuming, so we do it once)
27 1
    net["_is_elements"] = aux._select_is_elements_numba(net, sequence=sequence)
28

29 1
    ppc = _init_ppc(net, sequence)
30

31 1
    _build_bus_ppc(net, ppc)
32 1
    _build_gen_ppc(net, ppc)
33 1
    _add_ext_grid_sc_impedance_zero(net, ppc)
34 1
    _build_branch_ppc_zero(net, ppc)
35

36
    # adds auxilary buses for open switches at branches
37 1
    _switch_branches(net, ppc)
38

39
    # add auxilary buses for out of service buses at in service lines.
40
    # Also sets lines out of service if they are connected to two out of service buses
41 1
    _branches_with_oos_buses(net, ppc)
42 1
    if hasattr(net, "_isolated_buses"):
43 1
        ppc["bus"][net._isolated_buses, 1] = 4.
44
    # generates "internal" ppci format (for powerflow calc) from "external" ppc format and updates the bus lookup
45
    # Note: Also reorders buses and gens in ppc
46 1
    ppci = _ppc2ppci(ppc, net)
47
    #net._ppc0 = ppc    <--Obsolete. now covered in _init_ppc
48 1
    return ppc, ppci
49

50 1
def _build_branch_ppc_zero(net, ppc):
51
    """
52
    Takes the empty ppc network and fills it with the zero imepdance branch values. The branch
53
    datatype will be np.complex 128 afterwards.
54

55
    .. note:: The order of branches in the ppc is:
56
            1. Lines
57
            2. Transformers
58

59
    **INPUT**:
60
        **net** -The pandapower format network
61

62
        **ppc** - The PYPOWER format network to fill in values
63

64
    """
65 1
    length = _initialize_branch_lookup(net)
66 1
    lookup = net._pd2ppc_lookups["branch"]
67 1
    mode = net._options["mode"]
68 1
    ppc["branch"] = np.zeros(shape=(length, branch_cols), dtype=np.complex128)
69 1
    if mode == "sc":
70 1
        from pandapower.shortcircuit.idx_brch import branch_cols_sc
71 1
        branch_sc = np.empty(shape=(length, branch_cols_sc), dtype=float)
72 1
        branch_sc.fill(np.nan)
73 1
        ppc["branch"] = np.hstack((ppc["branch"], branch_sc ))
74 1
    ppc["branch"][:, :13] = np.array([0, 0, 0, 0, 0, 250, 250, 250, 1, 0, 1, -360, 360])
75 1
    _add_line_sc_impedance_zero(net, ppc)
76 1
    _add_trafo_sc_impedance_zero(net, ppc)
77 1
    if "trafo3w" in lookup:
78 0
        raise NotImplementedError("Three winding transformers are not implemented for unbalanced calculations")
79

80

81 1
def _add_trafo_sc_impedance_zero(net, ppc, trafo_df=None):
82 1
    if trafo_df is None:
83 1
        trafo_df = net["trafo"]
84 1
    branch_lookup = net["_pd2ppc_lookups"]["branch"]
85 1
    if not "trafo" in branch_lookup:
86 1
        return
87 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
88 1
    mode = net["_options"]["mode"]
89 1
    trafo_model = net["_options"]["trafo_model"]
90 1
    f, t = branch_lookup["trafo"]
91 1
    trafo_df["_ppc_idx"] = range(f, t)
92 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
93 1
    buses_all, gs_all, bs_all = np.array([], dtype=int), np.array([]), \
94
                                np.array([])
95 1
    if not "vector_group" in trafo_df:
96 0
        raise ValueError("Vector Group of transformer needs to be specified for zero \
97
                         sequence modelling \n Try : net.trafo[\"vector_group\"] = 'Dyn'" )
98

99 1
    for vector_group, trafos in trafo_df.groupby("vector_group"):
100 1
        ppc_idx = trafos["_ppc_idx"].values.astype(int)
101 1
        ppc["branch"][ppc_idx, BR_STATUS] = 0
102

103 1
        if vector_group in ["Yy", "Yd", "Dy", "Dd"]:
104 1
            continue
105

106 1
        vk_percent = trafos["vk_percent"].values.astype(float)
107 1
        vkr_percent = trafos["vkr_percent"].values.astype(float)
108 1
        sn_mva = trafos["sn_mva"].values.astype(float)
109
        # Just put pos seq parameter if zero seq parameter is zero
110 1
        if not "vk0_percent" in trafos:
111 0
            raise ValueError("Short circuit voltage of transformer Vk0 needs to be specified for zero \
112
                             sequence modelling \n Try : net.trafo[\"vk0_percent\"] = net.trafo[\"vk_percent\"]" )
113 1
        vk0_percent = trafos["vk0_percent"].values.astype(float) if \
114
            trafos["vk0_percent"].values.astype(float).all() != 0. else \
115
            trafos["vk_percent"].values.astype(float)
116
        # Just put pos seq parameter if zero seq parameter is zero
117 1
        if not "vkr0_percent" in trafos:
118 0
            raise ValueError("Real part of short circuit voltage Vk0(Real) needs to be specified for transformer \
119
                             modelling \n Try : net.trafo[\"vkr0_percent\"] = net.trafo[\"vkr_percent\"]" )        
120 1
        vkr0_percent = trafos["vkr0_percent"].values.astype(float) if \
121
            trafos["vkr0_percent"].values.astype(float).all() != 0. else \
122
            trafos["vkr_percent"].values.astype(float)
123 1
        lv_buses = trafos["lv_bus"].values.astype(int)
124 1
        hv_buses = trafos["hv_bus"].values.astype(int)
125 1
        lv_buses_ppc = bus_lookup[lv_buses]
126 1
        hv_buses_ppc = bus_lookup[hv_buses]
127 1
        if not "mag0_percent" in trafos:
128
            # For Shell Type transformers vk0 = vk * 1
129
            #        and mag0_percent = 10 ... 100  Zm0/ Zsc0
130
            # --pg 50 DigSilent Power Factory Transformer manual
131 0
            raise ValueError("Magnetizing impedance to vk0 ratio needs to be specified for transformer \
132
                             modelling  \n Try : net.trafo[\"mag0_percent\"] = 100" )           
133 1
        mag0_ratio = trafos.mag0_percent.values.astype(float)
134 1
        if not "mag0_rx" in trafos:
135 0
            raise ValueError("Magnetizing impedance R/X ratio needs to be specified for transformer \
136
                             modelling \n Try : net.trafo[\"mag0_rx\"] = 0 " )        
137 1
        mag0_rx = trafos["mag0_rx"].values.astype(float)
138 1
        if not "si0_hv_partial" in trafos:
139 0
            raise ValueError("Zero sequence short circuit impedance partition towards HV side needs to be specified for transformer \
140
                             modelling \n Try : net.trafo[\"si0_hv_partial\"] = 0.9 " )          
141 1
        si0_hv_partial = trafos.si0_hv_partial.values.astype(float)
142 1
        parallel = trafos.parallel.values.astype(float)
143 1
        in_service = trafos["in_service"].astype(int)
144

145 1
        ppc["branch"][ppc_idx, F_BUS] = hv_buses_ppc
146 1
        ppc["branch"][ppc_idx, T_BUS] = lv_buses_ppc
147

148 1
        vn_trafo_hv, vn_trafo_lv, shift = _calc_tap_from_dataframe(net, trafos)
149 1
        vn_lv = ppc["bus"][lv_buses_ppc, BASE_KV]
150 1
        ratio = _calc_nominal_ratio_from_dataframe(ppc, trafos, vn_trafo_hv, \
151
                                                   vn_trafo_lv, bus_lookup)
152 1
        ppc["branch"][ppc_idx, TAP] = ratio
153 1
        ppc["branch"][ppc_idx, SHIFT] = shift
154

155
        # zero seq. transformer impedance
156 1
        tap_lv = np.square(vn_trafo_lv / vn_lv) * net.sn_mva
157 1
        if mode == 'pf_3ph':
158
            # =============================================================================
159
            #     Changing base from transformer base to Network base to get Zpu(Net)
160
            #     Zbase = (kV).squared/S_mva
161
            #     Zpu(Net)={Zpu(trafo) * Zb(trafo)} / {Zb(Net)}
162
            #        Note:
163
            #             Network base voltage is Line-Neutral voltage in each phase
164
            #             Line-Neutral voltage= Line-Line Voltage(vn_lv) divided by sq.root(3)
165
            # =============================================================================
166 1
            tap_lv = np.square(vn_trafo_lv / vn_lv) * (3 * net.sn_mva)
167

168 1
        z_sc = vk0_percent / 100. / sn_mva * tap_lv
169 1
        r_sc = vkr0_percent / 100. / sn_mva * tap_lv
170 1
        z_sc = z_sc.astype(float)
171 1
        r_sc = r_sc.astype(float)
172 1
        x_sc = np.sign(z_sc) * np.sqrt(z_sc ** 2 - r_sc ** 2)
173 1
        z0_k = (r_sc + x_sc * 1j) / parallel
174 1
        y0_k = 1 / z0_k #adding admittance for "pi" model
175 1
        if mode == "sc":# or trafo_model == "pi":
176 1
            from pandapower.shortcircuit.idx_bus import C_MAX
177 1
            cmax = net._ppc["bus"][lv_buses_ppc, C_MAX]
178 1
            kt = _transformer_correction_factor(vk_percent, vkr_percent, \
179
                                                sn_mva, cmax)
180 1
            z0_k *= kt
181 1
            y0_k = 1 / z0_k
182
        # =============================================================================
183
        #       Transformer magnetising impedance for zero sequence
184
        # =============================================================================
185 1
        z_m = z_sc * mag0_ratio
186 1
        x_m = z_m / np.sqrt(mag0_rx ** 2 + 1)
187 1
        r_m = x_m * mag0_rx
188 1
        r0_trafo_mag = r_m / parallel
189 1
        x0_trafo_mag = x_m / parallel
190 1
        z0_mag = r0_trafo_mag + x0_trafo_mag * 1j
191
        # =============================================================================
192
        #         Star - Delta conversion ( T model to Pi Model)
193
        #      ----------- |__zc=ZAB__|-----------------
194
        #            _|                   _|
195
        #     za=ZAN|_|                  |_| zb=ZBN
196
        #            |                    |
197
        # =============================================================================
198 1
        z1 = si0_hv_partial * z0_k
199 1
        z2 = (1 - si0_hv_partial) * z0_k
200 1
        z3 = z0_mag
201 1
        z_temp = z1 * z2 + z2 * z3 + z1 * z3 
202 1
        za = z_temp / z2
203
#        za = z_temp / (z2+z3)
204 1
        zb = z_temp / z1
205
#        zb = z_temp / (z1+z3)
206 1
        zc = z_temp / z3  # ZAB  Transfer impedance
207
#        zc = z_temp / (z1+z2)  # ZAB  Transfer impedance
208 1
        YAB = 1 / zc.astype(complex)
209 1
        YAN = 1 / za.astype(complex)
210 1
        YBN = 1 / zb.astype(complex)
211
        
212
#        YAB_AN = (zc + za) /(zc * za).astype(complex)  # Series conn YAB and YAN
213
#        YAB_BN = (zc + zb) / (zc * zb).astype(complex)  # Series conn YAB and YBN
214

215 1
        YAB_AN = 1 / (zc + za).astype(complex)  # Series conn YAB and YAN
216 1
        YAB_BN = 1 / (zc + zb).astype(complex)  # Series conn YAB and YBN
217

218 1
        if vector_group == "Dyn":
219 1
            buses_all = np.hstack([buses_all, lv_buses_ppc])
220 1
            if trafo_model == "pi":
221 1
                y = y0_k # pi model
222
            else:
223 1
                y = (YAB + YBN).astype(complex) * int(ppc["baseMVA"])  # T model
224 1
            gs_all = np.hstack([gs_all, y.real * in_service])
225 1
            bs_all = np.hstack([bs_all, y.imag * in_service])
226

227 1
        elif vector_group == "YNd":
228 1
            buses_all = np.hstack([buses_all, hv_buses_ppc])
229 1
            if trafo_model == "pi":
230 1
                y = y0_k # pi model
231
            else:
232 0
                y = (YAB_BN + YAN).astype(complex) * int(ppc["baseMVA"]) #T model
233 1
            gs_all = np.hstack([gs_all, y.real * in_service])
234 1
            bs_all = np.hstack([bs_all, y.imag * in_service])
235

236 1
        elif vector_group == "Yyn":
237 1
            buses_all = np.hstack([buses_all, lv_buses_ppc])
238 1
            if trafo_model == "pi":
239 1
                y = 1/(z0_mag+z0_k).astype(complex)* int(ppc["baseMVA"]) #pi model
240
            else:
241
#                y = (YAB_AN + YBN).astype(complex) * int(ppc["baseMVA"]) #T model
242 0
                y = (YAB + YAB_BN + YBN).astype(complex)* int(ppc["baseMVA"])  # T model
243

244 1
            gs_all = np.hstack([gs_all, y.real * in_service])
245 1
            bs_all = np.hstack([bs_all, y.imag * in_service])
246
            
247

248 1
        elif vector_group == "YNyn":
249 1
            ppc["branch"][ppc_idx, BR_STATUS] = in_service
250
            # zc = ZAB
251 1
            ppc["branch"][ppc_idx, BR_R] = zc.real
252 1
            ppc["branch"][ppc_idx, BR_X] = zc.imag
253

254 1
            buses_all = np.hstack([buses_all, hv_buses_ppc])
255 1
            gs_all = np.hstack([gs_all, YAN.real * in_service \
256
                                * int(ppc["baseMVA"])])
257 1
            bs_all = np.hstack([bs_all, YAN.imag * in_service \
258
                                * int(ppc["baseMVA"])])
259

260 1
            buses_all = np.hstack([buses_all, lv_buses_ppc])
261 1
            gs_all = np.hstack([gs_all, YBN.real * in_service \
262
                                * int(ppc["baseMVA"])])
263 1
            bs_all = np.hstack([bs_all, YBN.imag * in_service \
264
                                * int(ppc["baseMVA"])])
265

266 1
        elif vector_group == "YNy":
267 1
            buses_all = np.hstack([buses_all, hv_buses_ppc])
268 1
            if trafo_model == "pi":
269 1
                y = 1/(z0_mag+z0_k).astype(complex)* int(ppc["baseMVA"])#pi model
270
            else:
271 0
                y = (YAB_BN + YAN).astype(complex) * int(ppc["baseMVA"])  #T model
272 1
            gs_all = np.hstack([gs_all, y.real * in_service])
273 1
            bs_all = np.hstack([bs_all, y.imag * in_service])
274

275 1
        elif vector_group == "Yzn":
276 1
            buses_all = np.hstack([buses_all, lv_buses_ppc])
277
            #            y = 1/(z0_mag+z0_k).astype(complex)* int(ppc["baseMVA"])#T model
278
            #            y= (za+zb+zc)/((za+zc)*zb).astype(complex)* int(ppc["baseMVA"])#pi model
279 1
            y = (YAB_AN + YBN).astype(complex) * int(ppc["baseMVA"])  #T model
280 1
            gs_all = np.hstack([gs_all, (1.1547) * y.real * in_service \
281
                                * int(ppc["baseMVA"])])
282 1
            bs_all = np.hstack([bs_all, (1.1547) * y.imag * in_service \
283
                                * int(ppc["baseMVA"])])
284

285 0
        elif vector_group[-1].isdigit():
286 0
            raise ValueError("Unknown transformer vector group %s -\
287
                             please specify vector group without \
288
                             phase shift number. Phase shift can be \
289
                             specified in net.trafo.shift_degree" % vector_group)
290
        else:
291 0
            raise ValueError("Transformer vector group %s is unknown\
292
                    / not implemented for three phase load flow" % vector_group)
293

294 1
    buses, gs, bs = aux._sum_by_group(buses_all, gs_all, bs_all)
295 1
    ppc["bus"][buses, GS] += gs
296 1
    ppc["bus"][buses, BS] += bs
297 1
    del net.trafo["_ppc_idx"]
298

299 1
def _add_ext_grid_sc_impedance_zero(net, ppc):
300 1
    mode = net["_options"]["mode"]
301

302 1
    if mode == "sc":
303 1
        from pandapower.shortcircuit.idx_bus import C_MAX, C_MIN
304 1
        case = net._options["case"]
305
    else:
306 1
        case = "max"
307 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
308 1
    eg = net["ext_grid"][net._is_elements["ext_grid"]]
309 1
    if len(eg) == 0:
310 0
        return
311 1
    eg_buses = eg.bus.values
312 1
    eg_buses_ppc = bus_lookup[eg_buses]
313

314 1
    if mode == "sc":
315 1
        c = ppc["bus"][eg_buses_ppc, C_MAX] if case == "max" else ppc["bus"][eg_buses_ppc, C_MIN]
316
    else:
317 1
        c = 1.1
318 1
    if not "s_sc_%s_mva" % case in eg:
319 0
        raise ValueError("short circuit apparent power s_sc_%s_mva needs to be specified for "% case +
320
                         "external grid" )
321 1
    s_sc = eg["s_sc_%s_mva" % case].values
322 1
    if not "rx_%s" % case in eg:
323 0
        raise ValueError("short circuit R/X rate rx_%s needs to be specified for external grid" %
324
                         case)
325 1
    rx = eg["rx_%s" % case].values
326

327 1
    z_grid = c / s_sc
328 1
    if mode == 'pf_3ph':
329 1
        z_grid = c / (s_sc/3)                       
330 1
    x_grid = z_grid / np.sqrt(rx ** 2 + 1)
331 1
    r_grid = rx * x_grid
332 1
    eg["r"] = r_grid
333 1
    eg["x"] = x_grid
334

335
    # ext_grid zero sequence impedance
336 1
    if case == "max":
337 1
        x0_grid = net.ext_grid["x0x_%s" % case] * x_grid
338 1
        r0_grid = net.ext_grid["r0x0_%s" % case] * x0_grid
339 1
    elif case == "min":
340 1
        x0_grid = net.ext_grid["x0x_%s" % case] * x_grid
341 1
        r0_grid = net.ext_grid["r0x0_%s" % case] * x0_grid
342 1
    y0_grid = 1 / (r0_grid + x0_grid*1j)
343

344 1
    buses, gs, bs = aux._sum_by_group(eg_buses_ppc, y0_grid.values.real, y0_grid.values.imag)
345 1
    ppc["bus"][buses, GS] = gs
346 1
    ppc["bus"][buses, BS] = bs
347

348

349 1
def _add_line_sc_impedance_zero(net, ppc):
350 1
    branch_lookup = net["_pd2ppc_lookups"]["branch"]
351 1
    mode = net["_options"]["mode"]                            
352 1
    if not "line" in branch_lookup:
353 0
        return
354 1
    line = net["line"]
355 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
356 1
    length = line["length_km"].values
357 1
    parallel = line["parallel"].values
358

359 1
    fb = bus_lookup[line["from_bus"].values]
360 1
    tb = bus_lookup[line["to_bus"].values]
361 1
    baseR = np.square(ppc["bus"][fb, BASE_KV]) / net.sn_mva
362 1
    if mode == 'pf_3ph':
363 1
        baseR = np.square(ppc["bus"][fb, BASE_KV]) / (3*net.sn_mva)                     
364 1
    f, t = branch_lookup["line"]
365
    # line zero sequence impedance
366 1
    ppc["branch"][f:t, F_BUS] = fb
367 1
    ppc["branch"][f:t, T_BUS] = tb
368 1
    ppc["branch"][f:t, BR_R] = line["r0_ohm_per_km"].values * length / baseR / parallel
369 1
    if mode == "sc":
370
        # temperature correction
371 1
        if net["_options"]["case"] == "min":
372 1
            ppc["branch"][f:t, BR_R] *= _end_temperature_correction_factor(net, short_circuit=True)
373 1
    ppc["branch"][f:t, BR_X] = line["x0_ohm_per_km"].values * length / baseR / parallel
374 1
    ppc["branch"][f:t, BR_B] = (2 * net["f_hz"] * math.pi * line["c0_nf_per_km"].values * 1e-9 * baseR * length * parallel)
375 1
    ppc["branch"][f:t, BR_STATUS] = line["in_service"].astype(int)

Read our documentation on viewing source code .

Loading