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

7 1
from collections import defaultdict, Counter
8 1
from itertools import chain
9

10 1
import numpy as np
11 1
import pandas as pd
12

13 1
from pandapower.auxiliary import _sum_by_group, phase_to_sequence
14 1
from pandapower.pypower.idx_bus import BUS_I, BASE_KV, PD, QD, GS, BS, VMAX, VMIN, BUS_TYPE, NONE, VM, VA, \
15
    CID, CZD, bus_cols, REF
16

17 1
try:
18 1
    from numba import jit
19 0
except ImportError:
20 0
    from .pf.no_numba import jit
21

22

23 1
@jit(nopython=True, cache=False)
24
def ds_find(ar, bus):  # pragma: no cover
25
    while True:
26
        p = ar[bus]
27
        if p == bus:
28
            break
29
        bus = p
30
    return p
31

32

33 1
@jit(nopython=True, cache=False)
34
def ds_union(ar, bus1, bus2, bus_is_pv):  # pragma: no cover
35
    root1 = ds_find(ar, bus1)
36
    root2 = ds_find(ar, bus2)
37
    if root1 == root2:
38
        return
39
    if bus_is_pv[root2]:
40
        ar[root1] = root2
41
    else:
42
        ar[root2] = root1
43

44

45 1
@jit(nopython=True, cache=False)
46
def ds_create(ar, switch_bus, switch_elm, switch_et_bus, switch_closed, switch_z_ohm,
47
              bus_is_pv, bus_in_service):  # pragma: no cover
48
    for i in range(len(switch_bus)):
49
        if not switch_closed[i] or not switch_et_bus[i] or switch_z_ohm[i] > 0:
50
            continue
51
        bus1 = switch_bus[i]
52
        bus2 = switch_elm[i]
53
        if bus_in_service[bus1] and bus_in_service[bus2]:
54
            ds_union(ar, bus1, bus2, bus_is_pv)
55

56

57 1
@jit(nopython=True, cache=False)
58 1
def fill_bus_lookup(ar, bus_lookup, bus_index):
59 0
    for i in range(len(bus_index)):
60 0
        bus_lookup[bus_index[i]] = i
61 0
    for b in bus_index:
62 0
        ds = ds_find(ar, b)
63 0
        bus_lookup[b] = bus_lookup[ar[ds]]
64

65

66 1
def create_bus_lookup_numba(net, bus_index, bus_is_idx, gen_is_mask, eg_is_mask):
67 1
    max_bus_idx = np.max(bus_index)
68
    # extract numpy arrays of switch table data
69 1
    switch = net["switch"]
70 1
    switch_bus = switch["bus"].values
71 1
    switch_elm = switch["element"].values
72 1
    switch_et_bus = switch["et"].values == "b"
73 1
    switch_closed = switch["closed"].values
74 1
    switch_z_ohm = switch['z_ohm'].values
75
    # create array for fast checking if a bus is in_service
76 1
    bus_in_service = np.zeros(max_bus_idx + 1, dtype=bool)
77 1
    bus_in_service[bus_is_idx] = True
78
    # create array for fast checking if a bus is pv bus
79 1
    bus_is_pv = np.zeros(max_bus_idx + 1, dtype=bool)
80 1
    bus_is_pv[net["ext_grid"]["bus"].values[eg_is_mask]] = True
81 1
    bus_is_pv[net["gen"]["bus"].values[gen_is_mask]] = True
82
    # create array that represents the disjoint set
83 1
    ar = np.arange(max_bus_idx + 1)
84 1
    ds_create(ar, switch_bus, switch_elm, switch_et_bus, switch_closed, switch_z_ohm, bus_is_pv, bus_in_service)
85
    # finally create and fill bus lookup
86 1
    bus_lookup = -np.ones(max_bus_idx + 1, dtype=int)
87 1
    fill_bus_lookup(ar, bus_lookup, bus_index)
88 1
    return bus_lookup
89

90

91 1
class DisjointSet(dict):
92 1
    def add(self, item):
93 0
        self[item] = item
94

95 1
    def find(self, item):
96 1
        parent = self[item]
97 1
        if self[parent] != parent:
98 1
            parent = self.find(parent)
99 1
            self[item] = parent
100 1
        return parent
101

102 1
    def union(self, item1, item2):
103 1
        p1 = self.find(item1)
104 1
        p2 = self.find(item2)
105 1
        self[p1] = p2
106

107

108 1
def create_consecutive_bus_lookup(net, bus_index):
109
    # create a mapping from arbitrary pp-index to a consecutive index starting at zero (ppc-index)
110 1
    consec_buses = np.arange(len(bus_index))
111
    # bus_lookup as dict:
112
    # bus_lookup = dict(zip(bus_index, consec_buses))
113
    # bus lookup as mask from pandapower -> pypower
114 1
    bus_lookup = -np.ones(max(bus_index) + 1, dtype=int)
115 1
    bus_lookup[bus_index] = consec_buses
116 1
    return bus_lookup
117

118

119 1
def create_bus_lookup_numpy(net, bus_index, bus_is_idx, gen_is_mask, eg_is_mask, closed_bb_switch_mask):
120 1
    bus_lookup = create_consecutive_bus_lookup(net, bus_index)
121 1
    net._fused_bb_switches = closed_bb_switch_mask & (net["switch"]["z_ohm"].values <= 0)
122 1
    if net._fused_bb_switches.any():
123
        # Note: this might seem a little odd - first constructing a pp to ppc mapping without
124
        # fused busses and then update the entries. The alternative (to construct the final
125
        # mapping at once) would require to determine how to fuse busses and which busses
126
        # are not part of any opened bus-bus switches first. It turns out, that the latter takes
127
        # quite some time in the average usecase, where #busses >> #bus-bus switches.
128

129
        # Find PV / Slack nodes -> their bus must be kept when fused with a PQ node
130 1
        pv_list = [net["ext_grid"]["bus"].values[eg_is_mask], net["gen"]["bus"].values[gen_is_mask]]
131 1
        pv_ref = np.unique(np.hstack(pv_list))
132
        # get the pp-indices of the buses which are connected to a switch to be fused
133 1
        fbus = net["switch"]["bus"].values[net._fused_bb_switches]
134 1
        tbus = net["switch"]["element"].values[net._fused_bb_switches]
135

136
        # create a mapping to map each bus to itself at frist ...
137 1
        ds = DisjointSet({e: e for e in chain(fbus, tbus)})
138

139
        # ... to follow each bus along a possible chain of switches to its final bus and update the
140
        # map
141 1
        for f, t in zip(fbus, tbus):
142 1
            ds.union(f, t)
143

144
        # now we can find out how to fuse each bus by looking up the final bus of the chain they
145
        # are connected to
146 1
        v = defaultdict(set)
147 1
        for a in ds:
148 1
            v[ds.find(a)].add(a)
149
        # build sets of buses which will be fused
150 1
        disjoint_sets = [e for e in v.values() if len(e) > 1]
151

152
        # check if PV buses need to be fused
153
        # if yes: the sets with PV buses must be found (which is slow)
154
        # if no: the check can be omitted
155 1
        if any(i in fbus or i in tbus for i in pv_ref):
156
            # for every disjoint set
157 1
            for dj in disjoint_sets:
158
                # check if pv buses are in the disjoint set dj
159 1
                pv_buses_in_set = set(pv_ref) & dj
160 1
                nr_pv_bus = len(pv_buses_in_set)
161 1
                if nr_pv_bus == 0:
162
                    # no pv buses. Use any bus in dj
163 1
                    map_to = bus_lookup[dj.pop()]
164
                else:
165
                    # one pv bus. Get bus from pv_buses_in_set
166 1
                    map_to = bus_lookup[pv_buses_in_set.pop()]
167 1
                for bus in dj:
168
                    # update lookup
169 1
                    bus_lookup[bus] = map_to
170
        else:
171
            # no PV buses in set
172 1
            for dj in disjoint_sets:
173
                # use any bus in set
174 1
                map_to = bus_lookup[dj.pop()]
175 1
                for bus in dj:
176
                    # update bus lookup
177 1
                    bus_lookup[bus] = map_to
178 1
    return bus_lookup
179

180

181 1
def create_bus_lookup(net, bus_index, bus_is_idx, gen_is_mask, eg_is_mask, numba):
182 1
    switches_with_pos_z_ohm = net["switch"]["z_ohm"].values > 0
183 1
    if switches_with_pos_z_ohm.any() or numba == False:
184
        # if there are any closed bus-bus switches find them
185 1
        closed_bb_switch_mask = ((net["switch"]["closed"].values == True) &
186
                                 (net["switch"]["et"].values == "b") &
187
                                 np.in1d(net["switch"]["bus"].values, bus_is_idx) &
188
                                 np.in1d(net["switch"]["element"].values, bus_is_idx))
189

190 1
    if switches_with_pos_z_ohm.any():
191 1
        net._impedance_bb_switches = closed_bb_switch_mask & switches_with_pos_z_ohm
192
    else:
193 1
        net._impedance_bb_switches = np.zeros(switches_with_pos_z_ohm.shape)
194

195 1
    if numba:
196 1
        bus_lookup = create_bus_lookup_numba(net, bus_index, bus_is_idx,
197
                                             gen_is_mask, eg_is_mask)
198
    else:
199 1
        bus_lookup = create_bus_lookup_numpy(net, bus_index, bus_is_idx,
200
                                             gen_is_mask, eg_is_mask, closed_bb_switch_mask)
201 1
    return bus_lookup
202

203

204 1
def get_voltage_init_vector(net, init_v, mode, sequence=None):
205 1
    if isinstance(init_v, str):
206 1
        if init_v == "results":
207 1
            res_table = "res_bus" if sequence is None else "res_bus_3ph"
208

209
            # init voltage possible if bus results are available
210 1
            if res_table not in net or not net[res_table].index.equals(net.bus.index):
211
                # cannot init from results, since sorting of results is different from element table
212
                # TO BE REVIEWED! Why there was no raise before this commit?
213 0
                raise UserWarning("Init from results not possible. Index of %s do not match with bus. "
214
                                  "You should sort res_bus before calling runpp." % res_table)
215

216 1
            if res_table == "res_bus_3ph":
217 1
                vm = net.res_bus_3ph[["vm_a_pu", "vm_b_pu", "vm_c_pu"]].values.T
218 1
                va = net.res_bus_3ph[["va_a_degree", "va_b_degree", "va_c_degree"]].values.T
219

220 1
                voltage_vector = phase_to_sequence(vm * np.exp(1j * np.pi * va / 180.))[sequence, :]
221

222 1
                if mode == "magnitude":
223 1
                    return np.abs(voltage_vector)
224 1
                elif mode == "angle":
225 1
                    return np.angle(voltage_vector) * 180 / np.pi
226
                else:
227 0
                    raise UserWarning(str(mode)+" for initialization not available!")
228
            else:
229 1
                if mode == "magnitude":
230 1
                    return net[res_table]["vm_pu"].values.copy()
231 1
                elif mode == "angle":
232 1
                    return net[res_table]["va_degree"].values.copy()
233
                else:
234 0
                    raise UserWarning(str(mode)+" for initialization not available!")
235 1
        if init_v == "flat":
236 1
            return None
237 1
    elif isinstance(init_v, (float, np.ndarray, list)) and sequence is None or sequence == 1:
238 1
        return init_v
239 1
    elif isinstance(init_v, pd.Series) and sequence is None or sequence == 1:
240 0
        if init_v.index.equals(net.bus.index):
241 0
            return init_v.loc[net.bus.index]
242
        else:
243 0
            raise UserWarning("Voltage starting vector indices do not match bus indices")
244

245

246 1
def _build_bus_ppc(net, ppc, sequence=None):
247
    """
248
    Generates the ppc["bus"] array and the lookup pandapower indices -> ppc indices
249
    """
250 1
    init_vm_pu = net["_options"]["init_vm_pu"]
251 1
    init_va_degree = net["_options"]["init_va_degree"]
252 1
    mode = net["_options"]["mode"]
253 1
    numba = net["_options"]["numba"] if "numba" in net["_options"] else False
254

255
    # get bus indices
256 1
    nr_xward = len(net.xward)
257 1
    nr_trafo3w = len(net.trafo3w)
258 1
    aux = dict()
259 1
    if nr_xward > 0 or nr_trafo3w > 0:
260 1
        bus_indices = [net["bus"].index.values, np.array([], dtype=np.int64)]
261 1
        max_idx = max(net["bus"].index) + 1
262 1
        if nr_xward > 0:
263 1
            aux_xward = np.arange(max_idx, max_idx + nr_xward, dtype=np.int64)
264 1
            aux["xward"] = aux_xward
265 1
            bus_indices.append(aux_xward)
266 1
        if nr_trafo3w:
267 1
            aux_trafo3w = np.arange(max_idx + nr_xward, max_idx + nr_xward + nr_trafo3w)
268 1
            aux["trafo3w"] = aux_trafo3w
269 1
            bus_indices.append(aux_trafo3w)
270 1
        bus_index = np.concatenate(bus_indices)
271
    else:
272 1
        bus_index = net["bus"].index.values
273
    # get in service elements
274

275 1
    if mode == "nx":
276 1
        bus_lookup = create_consecutive_bus_lookup(net, bus_index)
277
    else:
278 1
        _is_elements = net["_is_elements"]
279 1
        eg_is_mask = _is_elements['ext_grid']
280 1
        gen_is_mask = _is_elements['gen']
281 1
        bus_is_idx = _is_elements['bus_is_idx']
282 1
        bus_lookup = create_bus_lookup(net, bus_index, bus_is_idx,
283
                                       gen_is_mask, eg_is_mask, numba=numba)
284

285 1
    n_bus_ppc = len(bus_index)
286
    # init ppc with empty values
287 1
    ppc["bus"] = np.zeros(shape=(n_bus_ppc, bus_cols), dtype=float)
288 1
    ppc["bus"][:, :15] = np.array([0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 1, 2, 0, 0., 0.])  # changes of
289
    # voltage limits (2 and 0) must be considered in check_opf_data
290

291 1
    if sequence is not None and sequence != 1:
292 1
        ppc["bus"][:, VM] = 0.
293

294 1
    if mode == "sc":
295 1
        from pandapower.shortcircuit.idx_bus import bus_cols_sc
296 1
        bus_sc = np.empty(shape=(n_bus_ppc, bus_cols_sc), dtype=float)
297 1
        bus_sc.fill(np.nan)
298 1
        ppc["bus"] = np.hstack((ppc["bus"], bus_sc))
299

300
    # apply consecutive bus numbers
301 1
    ppc["bus"][:, BUS_I] = np.arange(n_bus_ppc)
302

303 1
    n_bus = len(net.bus.index)
304
    # init voltages from net
305 1
    ppc["bus"][:n_bus, BASE_KV] = net["bus"]["vn_kv"].values
306
    # set buses out of service (BUS_TYPE == 4)
307 1
    if nr_xward > 0 or nr_trafo3w > 0:
308 1
        in_service = np.concatenate([net["bus"]["in_service"].values,
309
                                     net["xward"]["in_service"].values,
310
                                     net["trafo3w"]["in_service"].values])
311
    else:
312 1
        in_service = net["bus"]["in_service"].values
313 1
    ppc["bus"][~in_service, BUS_TYPE] = NONE
314 1
    if mode != "nx":
315 1
        set_reference_buses(net, ppc, bus_lookup, mode)
316 1
    vm_pu = get_voltage_init_vector(net, init_vm_pu, "magnitude", sequence=sequence)
317 1
    if vm_pu is not None:
318 1
        ppc["bus"][:n_bus, VM] = vm_pu
319

320 1
    va_degree = get_voltage_init_vector(net, init_va_degree, "angle", sequence=sequence)
321 1
    if va_degree is not None:
322 1
        ppc["bus"][:n_bus, VA] = va_degree
323

324 1
    if mode == "sc":
325 1
        _add_c_to_ppc(net, ppc)
326

327 1
    if net._options["mode"] == "opf":
328 1
        if "max_vm_pu" in net.bus:
329 1
            ppc["bus"][:n_bus, VMAX] = net["bus"].max_vm_pu.values
330
        else:
331 1
            ppc["bus"][:n_bus, VMAX] = 2.  # changes of VMAX must be considered in check_opf_data
332 1
        if "min_vm_pu" in net.bus:
333 1
            ppc["bus"][:n_bus, VMIN] = net["bus"].min_vm_pu.values
334
        else:
335 1
            ppc["bus"][:n_bus, VMIN] = 0.  # changes of VMIN must be considered in check_opf_data
336

337 1
    if len(net.xward):
338 1
        _fill_auxiliary_buses(net, ppc, bus_lookup, "xward", "bus", aux)
339

340 1
    if len(net.trafo3w):
341 1
        _fill_auxiliary_buses(net, ppc, bus_lookup, "trafo3w", "hv_bus", aux)
342 1
    net["_pd2ppc_lookups"]["bus"] = bus_lookup
343 1
    net["_pd2ppc_lookups"]["aux"] = aux
344

345

346 1
def _fill_auxiliary_buses(net, ppc, bus_lookup, element, bus_column, aux):
347 1
    element_bus_idx = bus_lookup[net[element][bus_column].values]
348 1
    aux_idx = bus_lookup[aux[element]]
349 1
    ppc["bus"][aux_idx, BASE_KV] = ppc["bus"][element_bus_idx, BASE_KV]
350 1
    if net._options["mode"] == "opf":
351 1
        ppc["bus"][aux_idx, VMIN] = ppc["bus"][element_bus_idx, VMIN]
352 1
        ppc["bus"][aux_idx, VMAX] = ppc["bus"][element_bus_idx, VMAX]
353

354 1
    if net._options["init_vm_pu"] == "results":
355 1
        ppc["bus"][aux_idx, VM] = net["res_%s" % element]["vm_internal_pu"].values
356
    else:
357 1
        ppc["bus"][aux_idx, VM] = ppc["bus"][element_bus_idx, VM]
358 1
    if net._options["init_va_degree"] == "results":
359 1
        ppc["bus"][aux_idx, VA] = net["res_%s" % element]["va_internal_degree"].values
360
    else:
361 1
        ppc["bus"][aux_idx, VA] = ppc["bus"][element_bus_idx, VA]
362

363

364 1
def set_reference_buses(net, ppc, bus_lookup, mode):
365 1
    if mode == "nx":
366 0
        return
367 1
    eg_buses = bus_lookup[net.ext_grid.bus.values[net._is_elements["ext_grid"]]]
368 1
    ppc["bus"][eg_buses, BUS_TYPE] = REF
369 1
    if mode == "sc":
370 1
        gen_slacks = net._is_elements["gen"]  # generators are slacks for short-circuit calculation
371
    else:
372 1
        gen_slacks = net._is_elements["gen"] & net.gen["slack"].values
373 1
    if gen_slacks.any():
374 1
        slack_buses = net.gen["bus"].values[gen_slacks]
375 1
        ppc["bus"][bus_lookup[slack_buses], BUS_TYPE] = REF
376

377

378 1
def _calc_pq_elements_and_add_on_ppc(net, ppc, sequence= None):
379
    # init values
380 1
    b, p, q = np.array([], dtype=int), np.array([]), np.array([])
381

382 1
    _is_elements = net["_is_elements"]
383 1
    voltage_depend_loads = net["_options"]["voltage_depend_loads"]
384 1
    mode = net["_options"]["mode"]
385 1
    pq_elements = ["load", "motor", "sgen", "storage", "ward", "xward"]
386 1
    for element in pq_elements:
387 1
        tab = net[element]
388 1
        if len(tab):
389 1
            if element == "load" and voltage_depend_loads:
390 1
                cz = tab["const_z_percent"].values / 100.
391 1
                ci = tab["const_i_percent"].values / 100.
392 1
                if ((cz + ci) > 1).any():
393 0
                    raise ValueError("const_z_percent + const_i_percent need to be less or equal to " +
394
                                     "100%!")
395

396
                # cumulative sum of constant-current loads
397 1
                b_zip = tab["bus"].values
398 1
                load_counter = Counter(b_zip)
399

400 1
                bus_lookup = net["_pd2ppc_lookups"]["bus"]
401 1
                b_zip = bus_lookup[b_zip]
402 1
                load_counter = {bus_lookup[k]: v for k, v in load_counter.items()}
403 1
                b_zip, ci_sum, cz_sum = _sum_by_group(b_zip, ci, cz)
404

405 1
                for bus, no_loads in load_counter.items():
406 1
                    ci_sum[b_zip == bus] /= no_loads
407 1
                    cz_sum[b_zip == bus] /= no_loads
408

409 1
                ppc["bus"][b_zip, CID] = ci_sum
410 1
                ppc["bus"][b_zip, CZD] = cz_sum
411 1
            active = _is_elements[element]
412 1
            sign = -1 if element == "sgen" else 1
413 1
            if element == "motor":
414 1
                p_mw, q_mvar = _get_motor_pq(net)
415 1
                p = np.hstack([p, p_mw])
416 1
                q = np.hstack([q, q_mvar])
417 1
            elif element.endswith("ward"):
418 1
                p = np.hstack([p, tab["ps_mw"].values * active * sign])
419 1
                q = np.hstack([q, tab["qs_mvar"].values * active * sign])
420
            else:
421 1
                scaling = tab["scaling"].values
422 1
                p = np.hstack([p, tab["p_mw"].values * active * scaling * sign])
423 1
                q = np.hstack([q, tab["q_mvar"].values * active * scaling * sign])
424 1
            b = np.hstack([b, tab["bus"].values])
425

426 1
    for element in ["asymmetric_load", "asymmetric_sgen"]:
427 1
        if len(net[element]) > 0 and mode == "pf":
428 1
            p_mw, q_mvar = _get_symmetric_pq_of_unsymetric_element(net, element)
429 1
            sign = -1 if element.endswith("sgen") else 1
430 1
            p = np.hstack([p, p_mw * sign])
431 1
            q = np.hstack([q, q_mvar * sign ])
432 1
            b = np.hstack([b, net[element]["bus"].values])
433

434
    # sum up p & q of bus elements
435 1
    if b.size:
436 1
        bus_lookup = net["_pd2ppc_lookups"]["bus"]
437 1
        b = bus_lookup[b]
438 1
        b, vp, vq = _sum_by_group(b, p, q)
439 1
        ppc["bus"][b, PD] = vp
440 1
        ppc["bus"][b, QD] = vq
441

442

443 1
def _get_symmetric_pq_of_unsymetric_element(net, element):
444 1
    scale = net["_is_elements"][element] * net[element]["scaling"].values.T
445 1
    q_mvar = np.sum(net[element][["q_a_mvar", "q_b_mvar", "q_c_mvar"]].values, axis=1)
446 1
    p_mw = np.sum(net[element][["p_a_mw", "p_b_mw", "p_c_mw"]].values, axis=1)
447 1
    return p_mw*scale, q_mvar*scale
448

449 1
def _get_motor_pq(net):
450 1
    tab = net["motor"]
451 1
    active = net._is_elements["motor"]
452 1
    scale = tab["loading_percent"].values/100 *tab["scaling"].values*active
453

454 1
    efficiency = tab["efficiency_percent"].values
455 1
    p_mech = tab["pn_mech_mw"].values
456 1
    cos_phi = tab["cos_phi"].values
457

458 1
    p_mw = p_mech / efficiency * 100 * scale
459 1
    s_mvar = p_mw / cos_phi
460 1
    q_mvar = np.sqrt(s_mvar**2 - p_mw**2)
461 1
    return p_mw, q_mvar
462

463 1
def _calc_shunts_and_add_on_ppc(net, ppc):
464
    # init values
465 1
    b, p, q = np.array([], dtype=int), np.array([]), np.array([])
466 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
467
    # get in service elements
468 1
    _is_elements = net["_is_elements"]
469

470 1
    s = net["shunt"]
471 1
    if len(s) > 0:
472 1
        vl = _is_elements["shunt"]
473 1
        v_ratio = (ppc["bus"][bus_lookup[s["bus"].values], BASE_KV] / s["vn_kv"].values) ** 2
474 1
        q = np.hstack([q, s["q_mvar"].values * s["step"].values * v_ratio * vl])
475 1
        p = np.hstack([p, s["p_mw"].values * s["step"].values * v_ratio * vl])
476 1
        b = np.hstack([b, s["bus"].values])
477

478 1
    w = net["ward"]
479 1
    if len(w) > 0:
480 1
        vl = _is_elements["ward"]
481 1
        q = np.hstack([q, w["qz_mvar"].values * vl])
482 1
        p = np.hstack([p, w["pz_mw"].values * vl])
483 1
        b = np.hstack([b, w["bus"].values])
484

485 1
    xw = net["xward"]
486 1
    if len(xw) > 0:
487 1
        vl = _is_elements["xward"]
488 1
        q = np.hstack([q, xw["qz_mvar"].values * vl])
489 1
        p = np.hstack([p, xw["pz_mw"].values * vl])
490 1
        b = np.hstack([b, xw["bus"].values])
491

492 1
    loss_location = net._options["trafo3w_losses"].lower()
493 1
    trafo3w = net["trafo3w"]
494 1
    if loss_location == "star" and len(trafo3w) > 0:
495 1
        pfe_mw = trafo3w["pfe_kw"].values * 1e-3
496 1
        i0 = trafo3w["i0_percent"].values
497 1
        sn_mva = trafo3w["sn_hv_mva"].values
498

499 1
        q_mvar = (sn_mva * i0 / 100.) ** 2 - pfe_mw ** 2
500 1
        q_mvar[q_mvar < 0] = 0
501 1
        q_mvar = np.sqrt(q_mvar)
502

503 1
        vn_hv_trafo = trafo3w["vn_hv_kv"].values
504 1
        vn_hv_bus = ppc["bus"][bus_lookup[trafo3w.hv_bus.values], BASE_KV]
505 1
        v_ratio = (vn_hv_bus / vn_hv_trafo) ** 2
506

507 1
        q = np.hstack([q, q_mvar * v_ratio])
508 1
        p = np.hstack([p, pfe_mw * v_ratio])
509 1
        b = np.hstack([b, net._pd2ppc_lookups["aux"]["trafo3w"]])
510

511
    # if array is not empty
512 1
    if b.size:
513 1
        b = bus_lookup[b]
514 1
        b, vp, vq = _sum_by_group(b, p, q)
515

516 1
        ppc["bus"][b, GS] = vp
517 1
        ppc["bus"][b, BS] = -vq
518

519
# Short circuit relevant routines
520 1
def _add_ext_grid_sc_impedance(net, ppc):
521 1
    from pandapower.shortcircuit.idx_bus import C_MAX, C_MIN
522 1
    mode = net._options["mode"]
523 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
524 1
    if mode == "sc":
525 1
        case = net._options["case"]
526
    else:
527 1
        case = "max"
528 1
    eg = net["ext_grid"][net._is_elements["ext_grid"]]
529 1
    if len(eg) == 0:
530 1
        return
531 1
    eg_buses = eg.bus.values
532 1
    eg_buses_ppc = bus_lookup[eg_buses]
533

534 1
    if mode == "sc":
535 1
        c = ppc["bus"][eg_buses_ppc, C_MAX] if case == "max" else ppc["bus"][eg_buses_ppc, C_MIN]
536
    else:
537 1
        c = 1.1
538 1
    if not "s_sc_%s_mva" % case in eg:
539 0
        raise ValueError("short circuit apparent power s_sc_%s_mva needs to be specified for "% case +
540
                          "external grid \n Try: net.ext_grid['s_sc_max_mva'] = 1000" )
541 1
    s_sc = eg["s_sc_%s_mva" % case].values/ppc['baseMVA']
542 1
    if not "rx_%s" % case in eg:
543 0
        raise ValueError("short circuit R/X rate rx_%s needs to be specified for external grid \n Try: net.ext_grid['rx_max'] = 0.1" %
544
                          case)
545 1
    rx = eg["rx_%s" % case].values
546

547 1
    z_grid = c / s_sc
548 1
    if mode == 'pf_3ph':
549 1
        z_grid = c / (s_sc/3)  # 3 phase power divided to get 1 ph power
550 1
    x_grid = z_grid / np.sqrt(rx ** 2 + 1)
551 1
    r_grid = rx * x_grid
552 1
    eg["r"] = r_grid
553 1
    eg["x"] = x_grid
554

555 1
    y_grid = 1 / (r_grid + x_grid * 1j)
556 1
    buses, gs, bs = _sum_by_group(eg_buses_ppc, y_grid.real, y_grid.imag)
557 1
    ppc["bus"][buses, GS] = gs * ppc['baseMVA']
558 1
    ppc["bus"][buses, BS] = bs * ppc['baseMVA']
559 1
    return gs * ppc['baseMVA'], bs * ppc['baseMVA']
560

561 1
def _add_motor_impedances_ppc(net, ppc):
562 1
    if net._options["case"] == "min":
563 1
        return
564 1
    motor = net["motor"][net._is_elements["motor"]]
565 1
    if motor.empty:
566 1
        return
567 1
    for par in ["vn_kv", "lrc_pu", "efficiency_n_percent", "cos_phi_n", "rx", "pn_mech_mw"]:
568 1
        if any(pd.isnull(motor[par])):
569 0
            raise UserWarning("%s needs to be specified for all motors in net.motor.%s" % (par, par))
570 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
571 1
    motor_buses_ppc = bus_lookup[motor.bus.values]
572 1
    vn_net = ppc["bus"][motor_buses_ppc, BASE_KV]
573

574 1
    efficiency = motor.efficiency_n_percent.values
575 1
    cos_phi = motor.cos_phi_n.values
576 1
    p_mech = motor.pn_mech_mw.values
577 1
    vn_kv = motor.vn_kv.values
578 1
    lrc = motor.lrc_pu.values
579 1
    rx = motor.rx.values
580

581 1
    s_motor = p_mech / (efficiency/100 * cos_phi)
582 1
    z_motor_ohm = 1 / lrc * vn_kv**2 / s_motor
583 1
    z_motor_pu = z_motor_ohm / vn_net**2
584

585 1
    x_motor_pu = z_motor_pu / np.sqrt(rx ** 2 + 1)
586 1
    r_motor_pu = rx * x_motor_pu
587 1
    y_motor_pu = 1 / (r_motor_pu + x_motor_pu * 1j)
588

589 1
    buses, gs, bs = _sum_by_group(motor_buses_ppc, y_motor_pu.real, y_motor_pu.imag)
590 1
    ppc["bus"][buses, GS] = gs
591 1
    ppc["bus"][buses, BS] = bs
592

593

594 1
def _add_c_to_ppc(net, ppc):
595 1
    from pandapower.shortcircuit.idx_bus import C_MAX, C_MIN
596 1
    ppc["bus"][:, C_MAX] = 1.1
597 1
    ppc["bus"][:, C_MIN] = 1.
598 1
    lv_buses = np.where(ppc["bus"][:, BASE_KV] < 1.)
599 1
    if len(lv_buses) > 0:
600 1
        lv_tol_percent = net["_options"]["lv_tol_percent"]
601 1
        if lv_tol_percent == 10:
602 1
            c_ns = 1.1
603 1
        elif lv_tol_percent == 6:
604 1
            c_ns = 1.05
605
        else:
606 0
            raise ValueError("Voltage tolerance in the low voltage grid has" +
607
                             " to be either 6% or 10% according to IEC 60909")
608 1
        ppc["bus"][lv_buses, C_MAX] = c_ns
609 1
        ppc["bus"][lv_buses, C_MIN] = .95

Read our documentation on viewing source code .

Loading