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
# Use of this source code is governed by a BSD-style license that can be found in the LICENSE file.
6

7 1
import copy
8 1
import math
9 1
from functools import partial
10

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

14 1
from pandapower.auxiliary import get_values
15 1
from pandapower.pypower.idx_brch import F_BUS, T_BUS, BR_R, BR_X, BR_B, TAP, SHIFT, BR_STATUS, RATE_A, \
16
    BR_R_ASYM, BR_X_ASYM, branch_cols
17 1
from pandapower.pypower.idx_bus import BASE_KV, VM, VA
18

19

20 1
def _build_branch_ppc(net, ppc):
21
    """
22
    Takes the empty ppc network and fills it with the branch values. The branch
23
    datatype will be np.complex 128 afterwards.
24

25
    .. note:: The order of branches in the ppc is:
26
            1. Lines
27
            2. Transformers
28
            3. 3W Transformers (each 3W Transformer takes up three branches)
29
            4. Impedances
30
            5. Internal branch for extended ward
31

32
    **INPUT**:
33
        **net** -The pandapower format network
34

35
        **ppc** - The PYPOWER format network to fill in values
36

37
    """
38 1
    length = _initialize_branch_lookup(net)
39 1
    lookup = net._pd2ppc_lookups["branch"]
40 1
    mode = net._options["mode"]
41 1
    ppc["branch"] = np.zeros(shape=(length, branch_cols), dtype=np.complex128)
42
    # Check if this should be moved to somewhere else
43 1
    if mode == "sc":
44 1
        from pandapower.shortcircuit.idx_brch import branch_cols_sc
45 1
        branch_sc = np.empty(shape=(length, branch_cols_sc), dtype=float)
46 1
        branch_sc.fill(np.nan)
47 1
        ppc["branch"] = np.hstack((ppc["branch"], branch_sc))
48 1
    ppc["branch"][:, :13] = np.array([0, 0, 0, 0, 0, 250, 250, 250, 1, 0, 1, -360, 360])
49 1
    if "line" in lookup:
50 1
        _calc_line_parameter(net, ppc)
51 1
    if "trafo" in lookup:
52 1
        _calc_trafo_parameter(net, ppc)
53 1
    if "trafo3w" in lookup:
54 1
        _calc_trafo3w_parameter(net, ppc)
55 1
    if "impedance" in lookup:
56 1
        _calc_impedance_parameter(net, ppc)
57 1
    if "xward" in lookup:
58 1
        _calc_xward_parameter(net, ppc)
59 1
    if "switch" in lookup:
60 1
        _calc_switch_parameter(net, ppc)
61

62

63 1
def _initialize_branch_lookup(net):
64 1
    start = 0
65 1
    end = 0
66 1
    net._pd2ppc_lookups["branch"] = {}
67 1
    for element in ["line", "trafo", "trafo3w", "impedance", "xward"]:
68 1
        if len(net[element]) > 0:
69 1
            if element == "trafo3w":
70 1
                end = start + len(net[element]) * 3
71
            else:
72 1
                end = start + len(net[element])
73 1
            net._pd2ppc_lookups["branch"][element] = (start, end)
74 1
            start = end
75 1
    if "_impedance_bb_switches" in net and net._impedance_bb_switches.any():
76 1
        end = start + net._impedance_bb_switches.sum()
77 1
        net._pd2ppc_lookups["branch"]["switch"] = (start, end)
78 1
    return end
79

80

81 1
def _calc_trafo3w_parameter(net, ppc):
82 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
83 1
    branch = ppc["branch"]
84 1
    f, t = net["_pd2ppc_lookups"]["branch"]["trafo3w"]
85 1
    trafo_df = _trafo_df_from_trafo3w(net)
86 1
    hv_bus = get_trafo_values(trafo_df, "hv_bus").astype(int)
87 1
    lv_bus = get_trafo_values(trafo_df, "lv_bus").astype(int)
88 1
    in_service = get_trafo_values(trafo_df, "in_service").astype(int)
89 1
    branch[f:t, F_BUS] = bus_lookup[hv_bus]
90 1
    branch[f:t, T_BUS] = bus_lookup[lv_bus]
91 1
    r, x, y, ratio, shift = _calc_branch_values_from_trafo_df(net, ppc, trafo_df)
92 1
    branch[f:t, BR_R] = r
93 1
    branch[f:t, BR_X] = x
94 1
    branch[f:t, BR_B] = y
95 1
    branch[f:t, TAP] = ratio
96 1
    branch[f:t, SHIFT] = shift
97 1
    branch[f:t, BR_STATUS] = in_service
98 1
    if net["_options"]["mode"] == "opf":
99 1
        if "max_loading_percent" in trafo_df:
100 1
            max_load = get_trafo_values(trafo_df, "max_loading_percent")
101 1
            sn_mva = get_trafo_values(trafo_df, "sn_mva")
102 1
            branch[f:t, RATE_A] = max_load / 100. * sn_mva
103
        else:
104 1
            branch[f:t, RATE_A] = np.nan
105

106

107 1
def _calc_line_parameter(net, ppc, elm="line", ppc_elm="branch"):
108
    """
109
    calculates the line parameter in per unit.
110

111
    **INPUT**:
112
        **net** - The pandapower format network
113

114
        **ppc** - the ppc array
115

116
    **OPTIONAL**:
117
        **elm** - The pandapower element (normally "line")
118

119
        **ppc_elm** - The ppc element (normally "branch")
120

121
    **RETURN**:
122
        **t** - Temporary line parameter. Which is a complex128
123
                Nunmpy array. with the following order:
124
                0:bus_a; 1:bus_b; 2:r_pu; 3:x_pu; 4:b_pu
125
    """
126 1
    f, t = net._pd2ppc_lookups[ppc_elm][elm]
127 1
    branch = ppc[ppc_elm]
128 1
    mode = net["_options"]["mode"]
129 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
130 1
    line = net[elm]
131 1
    from_bus = bus_lookup[line["from_bus"].values]
132 1
    to_bus = bus_lookup[line["to_bus"].values]
133 1
    length_km = line["length_km"].values
134 1
    parallel = line["parallel"].values
135 1
    base_kv = ppc["bus"][from_bus, BASE_KV]
136 1
    baseR = np.square(base_kv) / (3 * net.sn_mva) if mode == "pf_3ph" else np.square(base_kv) / net.sn_mva 
137

138 1
    branch[f:t, F_BUS] = from_bus
139 1
    branch[f:t, T_BUS] = to_bus
140 1
    branch[f:t, BR_R] = line["r_ohm_per_km"].values * length_km / baseR / parallel
141 1
    branch[f:t, BR_X] = line["x_ohm_per_km"].values * length_km / baseR / parallel
142

143 1
    if mode == "sc":
144
        # temperature correction
145 1
        if net["_options"]["case"] == "min":
146 1
            branch[f:t, BR_R] *= _end_temperature_correction_factor(net, short_circuit=True)
147
    else:
148
        # temperature correction
149 1
        if net["_options"]["consider_line_temperature"]:
150 1
            branch[f:t, BR_R] *= _end_temperature_correction_factor(net)
151

152 1
        b = 2 * net.f_hz * math.pi * line["c_nf_per_km"].values * 1e-9 * baseR * length_km * parallel
153 1
        g = line["g_us_per_km"].values * 1e-6 * baseR * length_km * parallel
154 1
        branch[f:t, BR_B] = b - g * 1j
155
    # in service of lines
156 1
    branch[f:t, BR_STATUS] = line["in_service"].values
157 1
    if net._options["mode"] == "opf":
158
        # RATE_A is conisdered by the (PowerModels) OPF. If zero -> unlimited
159 1
        max_load = line.max_loading_percent.values if "max_loading_percent" in line else 0.
160 1
        vr = net.bus.loc[line["from_bus"].values, "vn_kv"].values * np.sqrt(3.)
161 1
        max_i_ka = line.max_i_ka.values
162 1
        df = line.df.values
163
        # This calculates the maximum apparent power at 1.0 p.u.
164 1
        branch[f:t, RATE_A] = max_load / 100. * max_i_ka * df * parallel * vr
165

166

167 1
def _calc_trafo_parameter(net, ppc):
168
    '''
169
    Calculates the transformer parameter in per unit.
170

171
    **INPUT**:
172
        **net** - The pandapower format network
173

174
    **RETURN**:
175
        **temp_para** -
176
        Temporary transformer parameter. Which is a np.complex128
177
        Numpy array. with the following order:
178
        0:hv_bus; 1:lv_bus; 2:r_pu; 3:x_pu; 4:b_pu; 5:tab, 6:shift
179
    '''
180

181 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
182 1
    f, t = net["_pd2ppc_lookups"]["branch"]["trafo"]
183 1
    branch = ppc["branch"]
184 1
    trafo = net["trafo"]
185 1
    parallel = trafo["parallel"].values
186 1
    branch[f:t, F_BUS] = bus_lookup[trafo["hv_bus"].values]
187 1
    branch[f:t, T_BUS] = bus_lookup[trafo["lv_bus"].values]
188 1
    r, x, y, ratio, shift = _calc_branch_values_from_trafo_df(net, ppc)
189 1
    branch[f:t, BR_R] = r
190 1
    branch[f:t, BR_X] = x
191 1
    branch[f:t, BR_B] = y
192 1
    branch[f:t, TAP] = ratio
193 1
    branch[f:t, SHIFT] = shift
194 1
    branch[f:t, BR_STATUS] = trafo["in_service"].values
195 1
    if any(trafo.df.values <= 0):
196 1
        raise UserWarning("Rating factor df must be positive. Transformers with false "
197
                          "rating factors: %s" % trafo.query('df<=0').index.tolist())
198 1
    if net._options["mode"] == "opf":
199 1
        max_load = trafo.max_loading_percent.values if "max_loading_percent" in trafo else 0
200 1
        sn_mva = trafo.sn_mva.values
201 1
        df = trafo.df.values
202 1
        branch[f:t, RATE_A] = max_load / 100. * sn_mva * df * parallel
203

204

205 1
def get_trafo_values(trafo_df, par):
206 1
    if isinstance(trafo_df, dict):
207 1
        return trafo_df[par]
208
    else:
209 1
        return trafo_df[par].values
210

211

212 1
def _calc_branch_values_from_trafo_df(net, ppc, trafo_df=None, sequence=1):
213
    """
214
    Calculates the MAT/PYPOWER-branch-attributes from the pandapower trafo dataframe.
215

216
    PYPOWER and MATPOWER uses the PI-model to model transformers.
217
    This function calculates the resistance r, reactance x, complex susceptance c and the tap ratio
218
    according to the given parameters.
219

220
    .. warning:: This function returns the subsceptance b as a complex number
221
        **(-img + -re*i)**. MAT/PYPOWER is only intended to calculate the
222
        imaginary part of the subceptance. However, internally c is
223
        multiplied by i. By using subsceptance in this way, it is possible
224
        to consider the ferromagnetic loss of the coil. Which would
225
        otherwise be neglected.
226

227

228
    .. warning:: Tab switches effect calculation as following:
229
        On **high-voltage** side(=1) -> only **tab** gets adapted.
230
        On **low-voltage** side(=2) -> **tab, x, r** get adapted.
231
        This is consistent with Sincal.
232
        The Sincal method in this case is questionable.
233

234

235
    **INPUT**:
236
        **pd_trafo** - The pandapower format Transformer Dataframe.
237
                        The Transformer modell will only readfrom pd_net
238

239
    **RETURN**:
240
        **temp_para** - Temporary transformer parameter. Which is a complex128
241
                        Nunmpy array. with the following order:
242
                        0:r_pu; 1:x_pu; 2:b_pu; 3:tab;
243

244
    """
245 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
246 1
    if trafo_df is None:
247 1
        trafo_df = net["trafo"]
248 1
    lv_bus = get_trafo_values(trafo_df, "lv_bus")
249 1
    vn_lv = ppc["bus"][bus_lookup[lv_bus], BASE_KV]
250
    ### Construct np.array to parse results in ###
251
    # 0:r_pu; 1:x_pu; 2:b_pu; 3:tab;
252 1
    vn_trafo_hv, vn_trafo_lv, shift = _calc_tap_from_dataframe(net, trafo_df)
253 1
    ratio = _calc_nominal_ratio_from_dataframe(ppc, trafo_df, vn_trafo_hv, vn_trafo_lv,
254
                                               bus_lookup)
255 1
    r, x, y = _calc_r_x_y_from_dataframe(net, trafo_df, vn_trafo_lv, vn_lv, ppc, sequence=sequence)
256 1
    return r, x, y, ratio, shift
257

258

259 1
def _calc_r_x_y_from_dataframe(net, trafo_df, vn_trafo_lv, vn_lv, ppc, sequence=1):
260 1
    mode = net["_options"]["mode"]
261 1
    trafo_model = net["_options"]["trafo_model"]
262

263 1
    r, x = _calc_r_x_from_dataframe(mode, trafo_df, vn_lv, vn_trafo_lv, net.sn_mva, sequence=sequence)
264

265 1
    if mode == "sc":
266 1
        y = 0
267 1
        if isinstance(trafo_df, pd.DataFrame):  # 2w trafo is dataframe, 3w trafo is dict
268 1
            from pandapower.shortcircuit.idx_bus import C_MAX
269 1
            bus_lookup = net._pd2ppc_lookups["bus"]
270 1
            cmax = ppc["bus"][bus_lookup[net.trafo.lv_bus.values], C_MAX]
271 1
            kt = _transformer_correction_factor(trafo_df.vk_percent, trafo_df.vkr_percent,
272
                                                trafo_df.sn_mva, cmax)
273 1
            r *= kt
274 1
            x *= kt
275
    else:
276 1
        y = _calc_y_from_dataframe(mode, trafo_df, vn_lv, vn_trafo_lv, net.sn_mva)
277

278 1
    if trafo_model == "pi":
279 1
        return r, x, y
280 1
    elif trafo_model == "t":
281 1
        return _wye_delta(r, x, y)
282
    else:
283 0
        raise ValueError("Unkonwn Transformer Model %s - valid values ar 'pi' or 't'" % trafo_model)
284

285

286 1
def _wye_delta(r, x, y):
287
    """
288
    20.05.2016 added by Lothar Löwer
289

290
    Calculate transformer Pi-Data based on T-Data
291

292
    """
293 1
    tidx = np.where(y != 0)
294 1
    za_star = (r[tidx] + x[tidx] * 1j) / 2
295 1
    zc_star = -1j / y[tidx]
296 1
    zSum_triangle = za_star * za_star + 2 * za_star * zc_star
297 1
    zab_triangle = zSum_triangle / zc_star
298 1
    zbc_triangle = zSum_triangle / za_star
299 1
    r[tidx] = zab_triangle.real
300 1
    x[tidx] = zab_triangle.imag
301 1
    y[tidx] = -2j / zbc_triangle
302 1
    return r, x, y
303

304

305 1
def _calc_y_from_dataframe(mode, trafo_df, vn_lv, vn_trafo_lv, sn_mva):
306
    """
307
    Calculate the subsceptance y from the transformer dataframe.
308

309
    INPUT:
310

311
        **trafo** (Dataframe) - The dataframe in net.trafo
312
        which contains transformer calculation values.
313

314
    OUTPUT:
315
        **subsceptance** (1d array, np.complex128) - The subsceptance in pu in
316
        the form (-b_img, -b_real)
317
    """
318

319 1
    baseR = np.square(vn_lv) / (3*sn_mva) if mode == 'pf_3ph' else np.square(vn_lv) / sn_mva
320 1
    vn_lv_kv = get_trafo_values(trafo_df, "vn_lv_kv")
321 1
    pfe = (get_trafo_values(trafo_df, "pfe_kw") * 1e-3) / 3 if mode == 'pf_3ph'\
322
        else get_trafo_values(trafo_df, "pfe_kw") * 1e-3
323 1
    parallel = get_trafo_values(trafo_df, "parallel")
324

325
    ### Calculate subsceptance ###
326 1
    vnl_squared = (vn_lv_kv ** 2)/3 if mode == 'pf_3ph' else vn_lv_kv ** 2
327 1
    b_real = pfe / vnl_squared * baseR
328 1
    i0 = get_trafo_values(trafo_df, "i0_percent") / 3 if mode == 'pf_3ph'\
329
        else get_trafo_values(trafo_df, "i0_percent")
330 1
    sn = get_trafo_values(trafo_df, "sn_mva")
331 1
    b_img = (i0 / 100. * sn) ** 2 - pfe ** 2
332

333 1
    b_img[b_img < 0] = 0
334 1
    b_img = np.sqrt(b_img) * baseR / vnl_squared
335 1
    y = - b_real * 1j - b_img * np.sign(i0)
336 1
    return y / np.square(vn_trafo_lv / vn_lv_kv) * parallel
337

338

339 1
def _calc_tap_from_dataframe(net, trafo_df):
340
    """
341
    Adjust the nominal voltage vnh and vnl to the active tab position "tap_pos".
342
    If "side" is 1 (high-voltage side) the high voltage vnh is adjusted.
343
    If "side" is 2 (low-voltage side) the low voltage vnl is adjusted
344

345
    INPUT:
346
        **net** - The pandapower format network
347

348
        **trafo** (Dataframe) - The dataframe in pd_net["structure"]["trafo"]
349
        which contains transformer calculation values.
350

351
    OUTPUT:
352
        **vn_hv_kv** (1d array, float) - The adusted high voltages
353

354
        **vn_lv_kv** (1d array, float) - The adjusted low voltages
355

356
        **trafo_shift** (1d array, float) - phase shift angle
357

358
    """
359 1
    calculate_voltage_angles = net["_options"]["calculate_voltage_angles"]
360 1
    mode = net["_options"]["mode"]
361 1
    vnh = copy.copy(get_trafo_values(trafo_df, "vn_hv_kv").astype(float))
362 1
    vnl = copy.copy(get_trafo_values(trafo_df, "vn_lv_kv").astype(float))
363 1
    trafo_shift = get_trafo_values(trafo_df, "shift_degree").astype(float) if calculate_voltage_angles else \
364
        np.zeros(len(vnh))
365 1
    if mode == "sc":
366 1
        return vnh, vnl, trafo_shift
367

368 1
    tap_pos = get_trafo_values(trafo_df, "tap_pos")
369 1
    tap_neutral = get_trafo_values(trafo_df, "tap_neutral")
370 1
    tap_diff = tap_pos - tap_neutral
371 1
    tap_phase_shifter = get_trafo_values(trafo_df, "tap_phase_shifter")
372 1
    tap_side = get_trafo_values(trafo_df, "tap_side")
373 1
    tap_step_percent = get_trafo_values(trafo_df, "tap_step_percent")
374 1
    tap_step_degree = get_trafo_values(trafo_df, "tap_step_degree")
375

376 1
    cos = lambda x: np.cos(np.deg2rad(x))
377 1
    sin = lambda x: np.sin(np.deg2rad(x))
378 1
    arctan = lambda x: np.rad2deg(np.arctan(x))
379

380 1
    for side, vn, direction in [("hv", vnh, 1), ("lv", vnl, -1)]:
381 1
        phase_shifters = tap_phase_shifter & (tap_side == side)
382 1
        tap_complex = np.isfinite(tap_step_percent) & np.isfinite(tap_pos) & (tap_side == side) & \
383
            ~phase_shifters
384 1
        if tap_complex.any():
385 1
            tap_steps = tap_step_percent[tap_complex] * tap_diff[tap_complex] / 100
386 1
            tap_angles = _replace_nan(tap_step_degree[tap_complex])
387 1
            u1 = vn[tap_complex]
388 1
            du = u1 * _replace_nan(tap_steps)
389 1
            vn[tap_complex] = np.sqrt((u1 + du * cos(tap_angles)) ** 2 + (du * sin(tap_angles)) ** 2)
390 1
            trafo_shift[tap_complex] += (arctan(direction * du * sin(tap_angles) /
391
                                                (u1 + du * cos(tap_angles))))
392 1
        if phase_shifters.any():
393 1
            degree_is_set = _replace_nan(tap_step_degree[phase_shifters]) != 0
394 1
            percent_is_set = _replace_nan(tap_step_percent[phase_shifters]) != 0
395 1
            if (degree_is_set & percent_is_set).any():
396 0
                raise UserWarning("Both tap_step_degree and tap_step_percent set for ideal phase shifter")
397 1
            trafo_shift[phase_shifters] += np.where(
398
                (degree_is_set),
399
                (direction * tap_diff[phase_shifters] * tap_step_degree[phase_shifters]),
400
                (direction * 2 * np.rad2deg(np.arcsin(tap_diff[phase_shifters] * \
401
                                                      tap_step_percent[phase_shifters] / 100 / 2)))
402
            )
403 1
    return vnh, vnl, trafo_shift
404

405

406 1
def _replace_nan(array, value=0):
407 1
    mask = np.isnan(array)
408 1
    array[mask] = value
409 1
    return array
410

411 1
def _calc_r_x_from_dataframe(mode, trafo_df, vn_lv, vn_trafo_lv, sn_mva, sequence=1):
412
    """
413
    Calculates (Vectorized) the resitance and reactance according to the
414
    transformer values
415

416
    """
417 1
    parallel = get_trafo_values(trafo_df, "parallel")
418 1
    if sequence == 1:
419 1
        vk_percent = get_trafo_values(trafo_df, "vk_percent")
420 1
        vkr_percent = get_trafo_values(trafo_df, "vkr_percent")
421 1
    elif sequence == 0:
422 1
        vk_percent = get_trafo_values(trafo_df, "vk0_percent")
423 1
        vkr_percent = get_trafo_values(trafo_df, "vkr0_percent")
424
    else:
425 0
        raise UserWarning("Unsupported sequence")
426 1
    tap_lv = np.square(vn_trafo_lv / vn_lv) * (3* sn_mva)  if mode == 'pf_3ph' else\
427
    np.square(vn_trafo_lv / vn_lv) * sn_mva  # adjust for low voltage side voltage converter
428

429 1
    sn_trafo_mva = get_trafo_values(trafo_df, "sn_mva")
430 1
    z_sc = vk_percent / 100. / sn_trafo_mva * tap_lv
431 1
    r_sc = vkr_percent / 100. / sn_trafo_mva * tap_lv
432 1
    x_sc = np.sign(z_sc) * np.sqrt((z_sc ** 2 - r_sc ** 2).astype(float))
433 1
    return r_sc / parallel, x_sc / parallel
434

435

436 1
def _calc_nominal_ratio_from_dataframe(ppc, trafo_df, vn_hv_kv, vn_lv_kv, bus_lookup):
437
    """
438
    Calculates (Vectorized) the off nominal tap ratio::
439

440
                  (vn_hv_kv / vn_lv_kv) / (ub1_in_kv / ub2_in_kv)
441

442
    INPUT:
443
        **net** (Dataframe) - The net for which to calc the tap ratio.
444

445
        **vn_hv_kv** (1d array, float) - The adjusted nominal high voltages
446

447
        **vn_lv_kv** (1d array, float) - The adjusted nominal low voltages
448

449
    OUTPUT:
450
        **tab** (1d array, float) - The off-nominal tap ratio
451
    """
452
    # Calculating tab (trasformer off nominal turns ratio)
453 1
    tap_rat = vn_hv_kv / vn_lv_kv
454 1
    hv_bus = get_trafo_values(trafo_df, "hv_bus")
455 1
    lv_bus = get_trafo_values(trafo_df, "lv_bus")
456 1
    nom_rat = get_values(ppc["bus"][:, BASE_KV], hv_bus, bus_lookup) / \
457
        get_values(ppc["bus"][:, BASE_KV], lv_bus, bus_lookup)
458 1
    return tap_rat / nom_rat
459

460

461 1
def _calc_impedance_parameter(net, ppc):
462 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
463 1
    f, t = net["_pd2ppc_lookups"]["branch"]["impedance"]
464 1
    branch = ppc["branch"]
465 1
    rij, xij, r_asym, x_asym = _calc_impedance_parameters_from_dataframe(net)
466 1
    branch[f:t, BR_R] = rij
467 1
    branch[f:t, BR_X] = xij
468 1
    branch[f:t, BR_R_ASYM] = r_asym
469 1
    branch[f:t, BR_X_ASYM] = x_asym
470 1
    branch[f:t, F_BUS] = bus_lookup[net.impedance["from_bus"].values]
471 1
    branch[f:t, T_BUS] = bus_lookup[net.impedance["to_bus"].values]
472 1
    branch[f:t, BR_STATUS] = net["impedance"]["in_service"].values
473

474

475 1
def _calc_impedance_parameters_from_dataframe(net):
476 1
    impedance = net.impedance
477 1
    sn_impedance = impedance["sn_mva"].values
478 1
    sn_net = net.sn_mva
479 1
    rij = impedance["rft_pu"].values
480 1
    xij = impedance["xft_pu"].values
481 1
    rji = impedance["rtf_pu"].values
482 1
    xji = impedance["xtf_pu"].values
483

484 1
    r = rij / sn_impedance * sn_net
485 1
    x = xij / sn_impedance * sn_net
486 1
    r_asym = (rji - rij) / sn_impedance * sn_net
487 1
    x_asym = (xji - xij) / sn_impedance * sn_net
488 1
    return r, x, r_asym, x_asym
489

490

491 1
def _calc_xward_parameter(net, ppc):
492 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
493 1
    f, t = net["_pd2ppc_lookups"]["branch"]["xward"]
494 1
    branch = ppc["branch"]
495 1
    baseR = np.square(get_values(ppc["bus"][:, BASE_KV], net["xward"]["bus"].values, bus_lookup)) / \
496
        net.sn_mva
497 1
    xw_is = net["_is_elements"]["xward"]
498 1
    branch[f:t, F_BUS] = bus_lookup[net["xward"]["bus"].values]
499 1
    branch[f:t, T_BUS] = bus_lookup[net._pd2ppc_lookups["aux"]["xward"]]
500 1
    branch[f:t, BR_R] = net["xward"]["r_ohm"] / baseR
501 1
    branch[f:t, BR_X] = net["xward"]["x_ohm"] / baseR
502 1
    branch[f:t, BR_STATUS] = xw_is
503

504

505 1
def _gather_branch_switch_info(bus, branch_id, branch_type, net):
506
    # determine at which end the switch is located
507
    # 1 = to-bus/lv-bus; 0 = from-bus/hv-bus
508 1
    branch_id = int(branch_id)
509 1
    lookup = net._pd2ppc_lookups["branch"]
510 1
    if branch_type == "l":
511 1
        side = "to" if net["line"]["to_bus"].at[branch_id] == bus else "from"
512 1
        branch_idx = net["line"].index.get_loc(branch_id)
513 1
        return side, int(bus), int(branch_idx)
514 1
    elif branch_type == "t":
515 1
        side = "hv" if net["trafo"]["hv_bus"].at[branch_id] == bus else "lv"
516 1
        branch_idx = lookup["trafo"][0] + net["trafo"].index.get_loc(branch_id)
517 1
        return side, int(bus), int(branch_idx)
518 1
    elif branch_type == "t3":
519 1
        f, t = lookup["trafo3w"]
520 1
        if net["trafo3w"]["hv_bus"].at[branch_id] == bus:
521 1
            side = "hv"
522 1
            offset = 0
523 1
        elif net["trafo3w"]["mv_bus"].at[branch_id] == bus:
524 1
            side = "mv"
525 1
            offset = (t - f) / 3
526 1
        elif net["trafo3w"]["lv_bus"].at[branch_id] == bus:
527 1
            side = "lv"
528 1
            offset = (t - f) / 3 * 2
529 1
        branch_idx = lookup["trafo3w"][0] + net["trafo3w"].index.get_loc(branch_id) + offset
530 1
        return side, int(bus), int(branch_idx)
531

532

533 1
def _switch_branches(net, ppc):
534 1
    from pandapower.shortcircuit.idx_bus import C_MIN, C_MAX
535 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
536 1
    calculate_voltage_angles = net._options["calculate_voltage_angles"]
537 1
    neglect_open_switch_branches = net._options["neglect_open_switch_branches"]
538 1
    mode = net._options["mode"]
539 1
    open_switches = (net.switch.closed.values == False)
540 1
    n_bus = ppc["bus"].shape[0]
541 1
    for et, element in [("l", "line"), ("t", "trafo"), ("t3", "trafo3w")]:
542 1
        switch_mask = open_switches & (net.switch.et.values == et)
543 1
        if not switch_mask.any():
544 1
            continue
545 1
        nr_open_switches = np.count_nonzero(switch_mask)
546 1
        mapfunc = partial(_gather_branch_switch_info, branch_type=et, net=net)
547 1
        switch_element = net["switch"]["element"].values[switch_mask]
548 1
        switch_buses = net["switch"]["bus"].values[switch_mask]
549 1
        switch_info = np.array(list(map(mapfunc, switch_buses, switch_element)))
550 1
        sw_sides = switch_info[:, 0]
551 1
        sw_bus_index = bus_lookup[switch_info[:, 1].astype(int)]
552 1
        sw_branch_index = switch_info[:, 2].astype(int)
553 1
        if neglect_open_switch_branches:
554
            # deactivate switches which have an open switch instead of creating aux buses
555 0
            ppc["branch"][sw_branch_index, BR_STATUS] = 0
556 0
            continue
557

558 1
        new_buses = np.zeros(shape=(nr_open_switches, ppc["bus"].shape[1]), dtype=float)
559 1
        new_buses[:, :15] = np.array([0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1.1, 0.9, 0, 0])
560 1
        new_indices = np.arange(n_bus, n_bus + nr_open_switches)
561 1
        new_buses[:, 0] = new_indices
562 1
        new_buses[:, BASE_KV] = ppc["bus"][sw_bus_index, BASE_KV]
563 1
        ppc["bus"] = np.vstack([ppc["bus"], new_buses])
564 1
        n_bus += new_buses.shape[0]
565 1
        init_vm = net._options["init_vm_pu"]
566 1
        init_va = net._options["init_va_degree"]
567 1
        for location in np.unique(sw_sides):
568 1
            mask = sw_sides == location
569 1
            buses = new_indices[mask]
570 1
            side = F_BUS if location == "hv" or location == "from" else T_BUS
571 1
            for init, col in [(init_vm, VM), (init_va, VA)]:
572 1
                if isinstance(init, str) and init == "results":
573 1
                    if col == VM:
574 1
                        res_column = net["res_%s" % element]["vm_%s_pu" % location]
575
                    else:
576 1
                        res_column = net["res_%s" % element]["va_%s_degree" % location]
577 1
                    init_values = res_column.loc[switch_element].values[mask]
578
                else:
579 1
                    if element == "line":
580 1
                        opposite_buses = ppc["branch"][sw_branch_index[mask], side].real.astype(int)
581 1
                        init_values = ppc["bus"][opposite_buses, col]
582
                    else:
583 1
                        opposite_side = T_BUS if side == F_BUS else F_BUS
584 1
                        opposite_buses = ppc["branch"][sw_branch_index[mask], opposite_side].real.astype(int)
585 1
                        if col == VM:
586 1
                            taps = ppc["branch"][sw_branch_index[mask], TAP].real
587 1
                            init_values = ppc["bus"][opposite_buses, col] * taps
588
                        else:
589 1
                            if calculate_voltage_angles:
590 1
                                shift = ppc["branch"][sw_branch_index[mask], SHIFT].real.astype(int)
591 1
                                init_values = ppc["bus"][opposite_buses, col] + shift
592
                            else:
593 1
                                init_values = ppc["bus"][opposite_buses, col]
594 1
                ppc["bus"][buses, col] = init_values
595 1
            if mode == "sc":
596 1
                ppc["bus"][buses, C_MAX] = ppc["bus"][opposite_buses, C_MAX]
597 1
                ppc["bus"][buses, C_MIN] = ppc["bus"][opposite_buses, C_MIN]
598 1
            ppc["branch"][sw_branch_index[mask], side] = new_indices[mask]
599

600

601 1
def _branches_with_oos_buses(net, ppc):
602
    """
603
    Updates the ppc["branch"] matrix with the changed from or to values
604
    if the branch is connected to an out of service bus
605

606
    Adds auxiliary buses if branch is connected to an out of service bus
607
    Sets branch out of service if connected to two out of service buses
608

609
    **INPUT**:
610
        **n** - The pandapower format network
611

612
        **ppc** - The PYPOWER format network to fill in values
613
        **bus_is** - The in service buses
614
    """
615 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
616
    # get in service elements
617 1
    _is_elements = net["_is_elements"]
618 1
    bus_is_idx = _is_elements['bus_is_idx']
619 1
    line_is_idx = _is_elements['line_is_idx']
620

621 1
    n_oos_buses = len(net['bus']) - len(bus_is_idx)
622

623
    # only filter lines at oos buses if oos buses exists
624 1
    if n_oos_buses > 0:
625 1
        n_bus = len(ppc["bus"])
626 1
        future_buses = [ppc["bus"]]
627
        # out of service buses
628 1
        bus_oos = np.setdiff1d(net['bus'].index.values, bus_is_idx)
629
        # from buses of line
630 1
        line_buses = net["line"][["from_bus", "to_bus"]].loc[line_is_idx].values
631 1
        f_bus = line_buses[:, 0]
632 1
        t_bus = line_buses[:, 1]
633

634
        # determine on which side of the line the oos bus is located
635 1
        mask_from = np.in1d(f_bus, bus_oos)
636 1
        mask_to = np.in1d(t_bus, bus_oos)
637

638 1
        mask_and = mask_to & mask_from
639 1
        if np.any(mask_and):
640 1
            mask_from[mask_and] = False
641 1
            mask_to[mask_and] = False
642

643
        # get lines that are connected to oos bus at exactly one side
644
        # buses that are connected to two oos buses will be removed by ext2int
645 1
        mask_or = mask_to | mask_from
646
        # check whether buses are connected to line
647 1
        oos_buses_at_lines = np.hstack([f_bus[mask_from], t_bus[mask_to]])
648 1
        n_oos_buses_at_lines = len(oos_buses_at_lines)
649

650
        # only if oos_buses are at lines (they could be isolated as well)
651 1
        if n_oos_buses_at_lines > 0:
652 1
            ls_info = np.zeros((n_oos_buses_at_lines, 3), dtype=int)
653 1
            ls_info[:, 0] = mask_to[mask_or] & ~mask_from[mask_or]
654 1
            ls_info[:, 1] = oos_buses_at_lines
655 1
            ls_info[:, 2] = np.nonzero(np.in1d(net['line'].index, line_is_idx[mask_or]))[0]
656

657
            # ls_info = list(map(mapfunc,
658
            #               line_switches["bus"].values,
659
            #               line_switches["element"].values))
660
            # we now have the following matrix
661
            # 0: 1 if switch is at to_bus, 0 else
662
            # 1: bus of the switch
663
            # 2: position of the line a switch is connected to
664
            # ls_info = np.array(ls_info, dtype=int)
665

666
            # build new buses
667 1
            new_ls_buses = np.zeros(shape=(n_oos_buses_at_lines, ppc["bus"].shape[1]), dtype=float)
668 1
            new_indices = np.arange(n_bus, n_bus + n_oos_buses_at_lines)
669
            # the newly created buses
670 1
            new_ls_buses[:, :15] = np.array([0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1.1, 0.9, 0, 0])
671 1
            new_ls_buses[:, 0] = new_indices
672 1
            new_ls_buses[:, BASE_KV] = get_values(ppc["bus"][:, BASE_KV], ls_info[:, 1], bus_lookup)
673

674 1
            future_buses.append(new_ls_buses)
675

676
            # re-route the end of lines to a new bus
677 1
            ppc["branch"][ls_info[ls_info[:, 0].astype(bool), 2], 1] = \
678
                new_indices[ls_info[:, 0].astype(bool)]
679 1
            ppc["branch"][ls_info[np.logical_not(ls_info[:, 0]), 2], 0] = \
680
                new_indices[np.logical_not(ls_info[:, 0])]
681

682 1
            ppc["bus"] = np.vstack(future_buses)
683

684

685 1
def _calc_switch_parameter(net, ppc):
686
    """
687
    calculates the line parameter in per unit.
688

689
    **INPUT**:
690
        **net** -The pandapower format network
691

692
    **RETURN**:
693
        **t** - Temporary line parameter. Which is a complex128
694
                Nunmpy array. with the following order:
695
                0:bus_a; 1:bus_b; 2:r_pu; 3:x_pu; 4:b_pu
696
    """
697 1
    rx_ratio = net["_options"]["switch_rx_ratio"]
698 1
    rz_ratio = rx_ratio / np.sqrt(1 + rx_ratio ** 2)
699 1
    xz_ratio = 1 / np.sqrt(1 + rx_ratio ** 2)
700

701 1
    f, t = net["_pd2ppc_lookups"]["branch"]["switch"]
702 1
    branch = ppc["branch"]
703 1
    bus_lookup = net["_pd2ppc_lookups"]["bus"]
704 1
    switch = net.switch[net._impedance_bb_switches]
705 1
    fb = bus_lookup[switch["bus"].values]
706 1
    tb = bus_lookup[switch["element"].values]
707 1
    baseR = np.square(ppc["bus"][fb, BASE_KV]) / net.sn_mva
708 1
    branch[f:t, F_BUS] = fb
709 1
    branch[f:t, T_BUS] = tb
710

711 1
    z_switch = switch['z_ohm'].values
712
    # x_switch will have the same value of r_switch to avoid zero dividence
713 1
    branch[f:t, BR_R] = z_switch / baseR * rz_ratio
714 1
    branch[f:t, BR_X] = z_switch / baseR * xz_ratio
715

716

717 1
def _end_temperature_correction_factor(net, short_circuit=False):
718
    """
719
    Function to calculate resistance correction factor for the given temperature ("endtemp_degree").
720
    When multiplied by the factor, the value of r_ohm_per_km will correspond to the resistance at
721
    the given temperature.
722

723
    In case of short circuit calculation, the relevant value for the temperature is
724
    "endtemp_degree", which stands for the final temperature of a line after the short circuit.
725
    The temperature coefficient "alpha" is a constant value of 0.004 in the short circuit
726
    calculation standard IEC 60909-0:2016.
727

728
    In case of a load flow calculation, the relelvant parameter is "temperature_degree_celsius",
729
    which is specified by the user and allows calculating load flow for a given operating
730
    temperature.
731

732
    The alpha value can be provided according to the used material for the purpose of load flow
733
    calculation, e.g. 0.0039 for copper or 0.00403 for aluminum. If alpha is not provided in the
734
    net.line table, the default value of 0.004 is used.
735

736
    The calculation of the electrical resistance is based on the formula R = R20(1+alpha*(T-20°C)),
737
    where R is the calculated resistance, R20 is the resistance at 20 °C, alpha is the temperature
738
    coefficient of resistance of the conducting material and T is the line temperature in °C.
739
    Accordingly, the resulting correction factor is (1+alpha*(T-20°C)).
740

741
    Args:
742
        net: pandapowerNet
743
        short_circuit: whether the factor is calculated in the scope of a short circuit calculation
744

745
    Returns:
746
        correction factor for line R, by which the line parameter should be multiplied to
747
                obtain the value of resistance at line temperature "endtemp_degree"
748

749
    """
750

751 1
    if short_circuit:
752
        # endtemp_degree is line temperature that is reached as the result of a short circuit
753
        # this value is the property of the lines
754 1
        if "endtemp_degree" not in net.line.columns:
755 0
            raise UserWarning("Specify end temperature for lines in net.line.endtemp_degree")
756

757 1
        delta_t_degree_celsius = net.line.endtemp_degree.values.astype(np.float64) - 20
758
        # alpha is the temperature correction factor for the electric resistance of the material
759
        # formula from standard, used this way in short-circuit calculation
760 1
        alpha = 4e-3
761
    else:
762
        # temperature_degree_celsius is line temperature for load flow calculation
763 1
        if "temperature_degree_celsius" not in net.line.columns:
764 0
            raise UserWarning("Specify line temperature in net.line.temperature_degree_celsius")
765

766 1
        delta_t_degree_celsius = net.line.temperature_degree_celsius.values.astype(np.float64) - 20
767

768 1
        if 'alpha' in net.line.columns:
769 1
            alpha = net.line.alpha.values.astype(np.float64)
770
        else:
771 1
            alpha = 4e-3
772

773 1
    r_correction_for_temperature = 1 + alpha * delta_t_degree_celsius
774

775 1
    return r_correction_for_temperature
776

777

778 1
def _transformer_correction_factor(vk, vkr, sn, cmax):
779
    """
780
        2W-Transformer impedance correction factor in short circuit calculations,
781
        based on the IEC 60909-0:2016 standard.
782
        Args:
783
            vk: transformer short-circuit voltage, percent
784
            vkr: real-part of transformer short-circuit voltage, percent
785
            sn: transformer rating, kVA
786
            cmax: voltage factor to account for maximum worst-case currents, based on the lv side
787

788
        Returns:
789
            kt: transformer impedance correction factor for short-circuit calculations
790
        """
791 1
    zt = vk / 100 / sn
792 1
    rt = vkr / 100 / sn
793 1
    xt = np.sqrt(zt ** 2 - rt ** 2)
794 1
    kt = 0.95 * cmax / (1 + .6 * xt * sn)
795 1
    return kt
796

797

798 1
def get_is_lines(net):
799
    """
800
    get indices of lines that are in service and save that information in net
801
    """
802 0
    _is_elements = net["_is_elements"]
803 0
    _is_elements["line"] = net["line"][net["line"]["in_service"].values.astype(bool)]
804

805

806 1
def _trafo_df_from_trafo3w(net, sequence=1):
807 1
    nr_trafos = len(net["trafo3w"])
808 1
    trafo2 = dict()
809 1
    sides = ["hv", "mv", "lv"]
810 1
    mode = net._options["mode"]
811 1
    loss_side = net._options["trafo3w_losses"].lower()
812 1
    nr_trafos = len(net["trafo3w"])
813 1
    t3 = net["trafo3w"]
814 1
    if sequence==1:
815 1
        _calculate_sc_voltages_of_equivalent_transformers(t3, trafo2, mode)
816 1
    elif sequence==0:
817 1
        if mode != "sc":
818 0
            raise NotImplementedError("0 seq impedance calculation only implemented for short-circuit calculation!")
819 1
        _calculate_sc_voltages_of_equivalent_transformers_zero_sequence(t3, trafo2,)
820
    else:
821 0
        raise UserWarning("Unsupported sequence for trafo3w convertion")
822 1
    _calculate_3w_tap_changers(t3, trafo2, sides)
823 1
    zeros = np.zeros(len(net.trafo3w))
824 1
    aux_buses = net._pd2ppc_lookups["aux"]["trafo3w"]
825 1
    trafo2["hv_bus"] = {"hv": t3.hv_bus.values, "mv": aux_buses, "lv": aux_buses}
826 1
    trafo2["lv_bus"] = {"hv": aux_buses, "mv": t3.mv_bus.values, "lv": t3.lv_bus.values}
827 1
    trafo2["in_service"] = {side: t3.in_service.values for side in sides}
828 1
    trafo2["i0_percent"] = {side: t3.i0_percent.values if loss_side == side else zeros for side in sides}
829 1
    trafo2["pfe_kw"] = {side: t3.pfe_kw.values if loss_side == side else zeros for side in sides}
830 1
    trafo2["vn_hv_kv"] = {side: t3.vn_hv_kv.values for side in sides}
831 1
    trafo2["vn_lv_kv"] = {side: t3["vn_%s_kv" % side].values for side in sides}
832 1
    trafo2["shift_degree"] = {"hv": np.zeros(nr_trafos), "mv": t3.shift_mv_degree.values,
833
                              "lv": t3.shift_lv_degree.values}
834 1
    trafo2["tap_phase_shifter"] = {side: np.zeros(nr_trafos).astype(bool) for side in sides}
835 1
    trafo2["parallel"] = {side: np.ones(nr_trafos) for side in sides}
836 1
    trafo2["df"] = {side: np.ones(nr_trafos) for side in sides}
837 1
    if net._options["mode"] == "opf" and "max_loading_percent" in net.trafo3w:
838 1
        trafo2["max_loading_percent"] = {side: net.trafo3w.max_loading_percent.values for side in sides}
839 1
    return {var: np.concatenate([trafo2[var][side] for side in sides]) for var in trafo2.keys()}
840

841

842 1
def _calculate_sc_voltages_of_equivalent_transformers(t3, t2, mode):
843 1
    vk_3w = np.stack([t3.vk_hv_percent.values, t3.vk_mv_percent.values, t3.vk_lv_percent.values])
844 1
    vkr_3w = np.stack([t3.vkr_hv_percent.values, t3.vkr_mv_percent.values, t3.vkr_lv_percent.values])
845 1
    sn = np.stack([t3.sn_hv_mva.values, t3.sn_mv_mva.values, t3.sn_lv_mva.values])
846

847 1
    vk_2w_delta = z_br_to_bus_vector(vk_3w, sn)
848 1
    vkr_2w_delta = z_br_to_bus_vector(vkr_3w, sn)
849 1
    if mode == "sc":
850 1
        kt = _transformer_correction_factor(vk_3w, vkr_3w, sn, 1.1)
851 1
        vk_2w_delta *= kt
852 1
        vkr_2w_delta *= kt
853 1
    vki_2w_delta = np.sqrt(vk_2w_delta ** 2 - vkr_2w_delta ** 2)
854 1
    vkr_2w = wye_delta_vector(vkr_2w_delta, sn)
855 1
    vki_2w = wye_delta_vector(vki_2w_delta, sn)
856 1
    vk_2w = np.sign(vki_2w) * np.sqrt(vki_2w ** 2 + vkr_2w ** 2)
857 1
    if np.any(vk_2w == 0):
858 0
        raise UserWarning("Equivalent transformer with zero impedance!")
859 1
    t2["vk_percent"] = {"hv": vk_2w[0, :], "mv": vk_2w[1, :], "lv": vk_2w[2, :]}
860 1
    t2["vkr_percent"] = {"hv": vkr_2w[0, :], "mv": vkr_2w[1, :], "lv": vkr_2w[2, :]}
861 1
    t2["sn_mva"] = {"hv": sn[0, :], "mv": sn[1, :], "lv": sn[2, :]}
862

863

864 1
def _calculate_sc_voltages_of_equivalent_transformers_zero_sequence(t3, t2):
865 1
    vk_3w = np.stack([t3.vk_hv_percent.values, t3.vk_mv_percent.values, t3.vk_lv_percent.values])
866 1
    vkr_3w = np.stack([t3.vkr_hv_percent.values, t3.vkr_mv_percent.values, t3.vkr_lv_percent.values])
867 1
    vk0_3w = np.stack([t3.vk0_hv_percent.values, t3.vk0_mv_percent.values, t3.vk0_lv_percent.values])
868 1
    vkr0_3w = np.stack([t3.vkr0_hv_percent.values, t3.vkr0_mv_percent.values, t3.vkr0_lv_percent.values])
869 1
    sn = np.stack([t3.sn_hv_mva.values, t3.sn_mv_mva.values, t3.sn_lv_mva.values])
870

871 1
    vk0_2w_delta = z_br_to_bus_vector(vk0_3w, sn)
872 1
    vkr0_2w_delta = z_br_to_bus_vector(vkr0_3w, sn)
873

874
    # Only for "sc", calculated with positive sequence value
875 1
    kt = _transformer_correction_factor(vk_3w, vkr_3w, sn, 1.1)
876 1
    vk0_2w_delta *= kt
877 1
    vkr0_2w_delta *= kt
878

879 1
    vki0_2w_delta = np.sqrt(vk0_2w_delta ** 2 - vkr0_2w_delta ** 2)
880 1
    vkr0_2w = wye_delta_vector(vkr0_2w_delta, sn)
881 1
    vki0_2w = wye_delta_vector(vki0_2w_delta, sn)
882 1
    vk0_2w = np.sign(vki0_2w) * np.sqrt(vki0_2w ** 2 + vkr0_2w ** 2)
883 1
    if np.any(vk0_2w == 0):
884 0
        raise UserWarning("Equivalent transformer with zero impedance!")
885 1
    t2["vk0_percent"] = {"hv": vk0_2w[0, :], "mv": vk0_2w[1, :], "lv": vk0_2w[2, :]}
886 1
    t2["vkr0_percent"] = {"hv": vkr0_2w[0, :], "mv": vkr0_2w[1, :], "lv": vkr0_2w[2, :]}
887 1
    t2["sn_mva"] = {"hv": sn[0, :], "mv": sn[1, :], "lv": sn[2, :]}
888

889

890 1
def z_br_to_bus_vector(z, sn):
891 1
    return sn[0, :] * np.array([z[0, :] / sn[[0, 1], :].min(axis=0), z[1, :] /
892
                                sn[[1, 2], :].min(axis=0), z[2, :] / sn[[0, 2], :].min(axis=0)])
893

894

895 1
def wye_delta(zbr_n, s):
896 0
    return .5 * s / s[0] * np.array([(zbr_n[0] + zbr_n[2] - zbr_n[1]),
897
                                     (zbr_n[1] + zbr_n[0] - zbr_n[2]),
898
                                     (zbr_n[2] + zbr_n[1] - zbr_n[0])])
899

900

901 1
def wye_delta_vector(zbr_n, s):
902 1
    return .5 * s / s[0, :] * np.array([(zbr_n[0, :] + zbr_n[2, :] - zbr_n[1, :]),
903
                                        (zbr_n[1, :] + zbr_n[0, :] - zbr_n[2, :]),
904
                                        (zbr_n[2, :] + zbr_n[1, :] - zbr_n[0, :])])
905

906

907 1
def _calculate_3w_tap_changers(t3, t2, sides):
908 1
    tap_variables = ["tap_side", "tap_pos", "tap_neutral", "tap_max", "tap_min", "tap_step_percent",
909
                     "tap_step_degree"]
910 1
    sides = ["hv", "mv", "lv"]
911 1
    nr_trafos = len(t3)
912 1
    empty = np.zeros(nr_trafos)
913 1
    empty.fill(np.nan)
914 1
    tap_arrays = {var: {side: empty.copy() for side in sides} for var in tap_variables}
915 1
    tap_arrays["tap_side"] = {side: np.array([None] * nr_trafos) for side in sides}
916 1
    at_star_point = t3.tap_at_star_point.values
917 1
    any_at_star_point = at_star_point.any()
918 1
    for side in sides:
919 1
        tap_mask = t3.tap_side.values == side
920 1
        for var in tap_variables:
921 1
            tap_arrays[var][side][tap_mask] = t3[var].values[tap_mask]
922

923
        # t3 trafos with tap changer at terminals
924 1
        tap_arrays["tap_side"][side][tap_mask] = "hv" if side == "hv" else "lv"
925

926
        # t3 trafos with tap changer at star points
927 1
        if any_at_star_point:
928 0
            mask_star_point = tap_mask & at_star_point
929 0
            tap_arrays["tap_side"][side][mask_star_point] = "lv" if side == "hv" else "hv"
930 0
            tap_arrays["tap_step_degree"][side][mask_star_point] += 180
931 1
    t2.update(tap_arrays)

Read our documentation on viewing source code .

Loading