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

3
# Copyright 1996-2015 PSERC. All rights reserved.
4
# Use of this source code is governed by a BSD-style
5
# license that can be found in the LICENSE file.
6

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

10

11 1
import numpy as np
12 1
from numba import jit
13 1
from scipy.sparse import csr_matrix, coo_matrix
14

15 1
from pandapower.pypower.idx_brch import F_BUS, T_BUS
16 1
from pandapower.pypower.idx_bus import GS, BS
17 1
from pandapower.pypower.makeYbus import branch_vectors
18

19

20 1
@jit(nopython=True, cache=False)
21
def gen_Ybus(Yf_x, Yt_x, Ysh, col_Y, f, t, f_sort, t_sort, nb, nl, r_nl):  # pragma: no cover
22
    """
23
    Fast calculation of Ybus
24
    """
25

26
    r_nb = range(nb)
27

28
    # allocate data of Ybus in CSR format
29
    # Note: More space is allocated than needed with empty.
30
    #       The matrix size will be reduced afterwards
31
    alloc_size = nl * 2 + nb
32
    Yx = np.empty(alloc_size, dtype=np.complex128)  # data
33
    Yp = np.zeros(nb + 1, dtype=np.int64)  # row pointer
34
    Yj = np.empty(alloc_size, dtype=np.int64)  # colum indices
35

36
    # index iterators
37
    # a = iterator of f, b = iterator of t, curRow = current Row
38
    a, b, curRow = 0, 0, 0
39
    # number of nonzeros (total), number of nonzeros per row
40
    nnz, nnz_row = 0, 0
41
    # flag checks if diagonal entry was added
42
    YshAdded = False
43

44
    for curRow in r_nb:
45
        nnz_row = 0
46
        # iterate rows of Ybus
47

48
        # add entries from Yf
49
        while a < nl and f[f_sort[a]] == curRow:
50
            # Entries from f_sort[a] in current row of Ybus
51
            for col in (r_nl[f_sort[a]], r_nl[f_sort[a]] + nl):
52
                # 'Has entry at column in Yf: %i ' % col
53
                if col_Y[col] == curRow and not YshAdded:
54
                    # add Ysh and Yf_x (diagonal element). If not already added
55
                    curVal = Yf_x[col] + Ysh[curRow]
56
                    YshAdded = True
57
                else:
58
                    # add only Yf_x
59
                    curVal = Yf_x[col]
60

61
                for k in range(Yp[curRow], Yp[curRow] + nnz_row):
62
                    if col_Y[col] == Yj[k]:
63
                        # if entry at column already exists add value
64
                        Yx[k] += curVal
65
                        break
66
                else:
67
                    # new entry in Ybus
68
                    Yx[nnz] = curVal
69
                    Yj[nnz] = col_Y[col]
70
                    nnz += 1
71
                    nnz_row += 1
72
            a += 1
73

74
        # add entries from Yt
75
        while b < nl and t[t_sort[b]] == curRow:
76
            # Entries from t_sort[b] in current row of Ybus
77
            for col in (r_nl[t_sort[b]], r_nl[t_sort[b]] + nl):
78
                # 'Has entry at column in Yt: %i ' % col
79
                if col_Y[col] == curRow and not YshAdded:
80
                    # add Ysh and Yf_x (diagonal element). If not already added
81
                    curVal = Yt_x[col] + Ysh[curRow]
82
                    YshAdded = True
83
                else:
84
                    # add only Yt_x
85
                    curVal = Yt_x[col]
86

87
                for k in range(Yp[curRow], Yp[curRow] + nnz_row):
88
                    if col_Y[col] == Yj[k]:
89
                        # if entry at column already exists add value
90
                        Yx[k] += curVal
91
                        break
92
                else:
93
                    # new entry in Ybus
94
                    Yx[nnz] = curVal
95
                    Yj[nnz] = col_Y[col]
96
                    nnz += 1
97
                    nnz_row += 1
98
            b += 1
99

100
        if not YshAdded:
101
            # check if diagonal entry was added. If not -> add if not zero
102
            if Ysh[curRow]:
103
                Yx[nnz] = Ysh[curRow]
104
                Yj[nnz] = curRow
105
                nnz += 1
106
                nnz_row += 1
107

108
        YshAdded = False
109
        # add number of nonzeros in row to row pointer
110
        Yp[curRow + 1] = nnz_row + Yp[curRow]
111
        curRow += 1
112

113
    return Yx, Yj, Yp, nnz
114

115

116 1
def makeYbus(baseMVA, bus, branch):
117
    """Builds the bus admittance matrix and branch admittance matrices.
118

119
    Returns the full bus admittance matrix (i.e. for all buses) and the
120
    matrices C{Yf} and C{Yt} which, when multiplied by a complex voltage
121
    vector, yield the vector currents injected into each line from the
122
    "from" and "to" buses respectively of each line. Does appropriate
123
    conversions to p.u.
124

125
    @see: L{makeSbus}
126

127
    @author: Ray Zimmerman (PSERC Cornell)
128
    @author: Richard Lincoln
129

130
    modified by Florian Schaefer (to use numba) (florian.schaefer@uni-kassel.de)
131
    """
132
    ## constants
133 1
    nb = bus.shape[0]  ## number of buses
134 1
    nl = branch.shape[0]  ## number of lines
135

136
    ## for each branch, compute the elements of the branch admittance matrix where
137
    ##
138
    ##      | If |   | Yff  Yft |   | Vf |
139
    ##      |    | = |          | * |    |
140
    ##      | It |   | Ytf  Ytt |   | Vt |
141
    ##
142 1
    Ytt, Yff, Yft, Ytf = branch_vectors(branch, nl)
143

144
    ## compute shunt admittance
145
    ## if Psh is the real power consumed by the shunt at V = 1.0 p.u.
146
    ## and Qsh is the reactive power injected by the shunt at V = 1.0 p.u.
147
    ## then Psh - j Qsh = V * conj(Ysh * V) = conj(Ysh) = Gs - j Bs,
148
    ## i.e. Ysh = Psh + j Qsh, so ...
149
    ## vector of shunt admittances
150 1
    Ysh = (bus[:, GS] + 1j * bus[:, BS]) / baseMVA
151

152
    ## build connection matrices
153 1
    f = np.real(branch[:, F_BUS]).astype(int)  ## list of "from" buses
154 1
    t = np.real(branch[:, T_BUS]).astype(int)  ## list of "to" buses
155

156
    ## build Yf and Yt such that Yf * V is the vector of complex branch currents injected
157
    ## at each branch's "from" bus, and Yt is the same for the "to" bus end
158 1
    i = np.hstack([np.arange(nl), np.arange(nl)])  ## double set of row indices
159

160 1
    Yf_x = np.hstack([Yff, Yft])
161 1
    Yt_x = np.hstack([Ytf, Ytt])
162 1
    col_Y = np.hstack([f, t])
163

164 1
    Yf = coo_matrix((Yf_x, (i, col_Y)), (nl, nb)).tocsr()
165 1
    Yt = coo_matrix((Yt_x, (i, col_Y)), (nl, nb)).tocsr()
166 1
    Yx, Yj, Yp, nnz = gen_Ybus(Yf_x, Yt_x, Ysh, col_Y, f, t, np.argsort(f), np.argsort(t), nb, nl,
167
                               np.arange(nl, dtype=np.int64))
168 1
    Ybus = csr_matrix((np.resize(Yx, nnz), np.resize(Yj, nnz), Yp), (nb, nb))
169 1
    return Ybus, Yf, Yt

Read our documentation on viewing source code .

Loading