Skip to content
GitLab
Commit 3db666e1 authored by Markus Mirz's avatar Markus Mirz
Browse files

restructuring powerflow

parent 78224a81
Hide whitespace changes
Inline Side-by-side
acs-state-estimation/powerflow.py → acs-state-estimation/bc_powerflow.py
View file @ 3db666e1
...@@ -4,9 +4,8 @@ import math ...@@ -4,9 +4,8 @@ import math
class Znod: class Znod:
def __init__(self, branches, nodes): def __init__(self, branches, nodes):
# Complex impedance matrix # Complex impedance matrix
self.Z = np.zeros( (nodes.num-1,branches.num), dtype=np.complex_ ) self.Z = np.zeros((nodes.num-1,branches.num), dtype=np.complex)
# Add impedance of each branch for k in range(1,nodes.num):
for k in range(1, nodes.num):
i = k-1 i = k-1
aa = branches.start[i]-2 aa = branches.start[i]-2
if aa == -1: if aa == -1:
...@@ -14,123 +13,112 @@ class Znod: ...@@ -14,123 +13,112 @@ class Znod:
else: else:
self.Z[i] = self.Z[aa] self.Z[i] = self.Z[aa]
self.Z[i][i] = branches.Z[i] self.Z[i][i] = branches.Z[i]
# Add one row of zeros # Add one row of zeros
zzz = np.zeros((1, branches.num)) zzz = np.zeros((1,branches.num))
self.Z = np.concatenate( (zzz,self.Z), axis=0 ) self.Z = np.concatenate((zzz,self.Z), axis=0)
# Resistance matrix
self.R = self.Z.real self.R = self.Z.real
# Reactance matrix
self.X = self.Z.imag self.X = self.Z.imag
def BC_power_flow(branch, node): def BC_power_flow(branches, nodes):
"""It performs Power Flow by using branch current state variables.""" """It performs Power Flow by using rectangular branches current state variables."""
znod = Znod(branch,node) znod = Znod(branches,nodes)
# real matrices storing complex numbers # real matrices storing complex numbers
# real in pos n and imag in pos n+1 # real in position n and imag in position n+1
z = np.zeros(2*(branch.num+1)) z = np.zeros(2*(branches.num+1))
h = np.zeros(2*(branch.num+1)) h = np.zeros(2*(branches.num+1))
H = np.zeros((2*(branch.num+1),2*(branch.num+1))) H = np.zeros((2*(branches.num+1),2*(branches.num+1)))
for k in range(1, node.num+1): for k in range(1,nodes.num+1):
# i is matrix index
# the size of z,h,H is twice the number of branches
# but here we are using the node number to work on the
# matrix. Is this not dangerous?
i = k-1 i = k-1
m = 2*i m = 2*i
if node.type[i] == 'slack': if nodes.type[i] == 'slack':
V = node.pwr_flow_values[1][i] V = nodes.pwr_flow_values[1][i]
theta = node.pwr_flow_values[2][i] theta = nodes.pwr_flow_values[2][i]
# real node voltage
z[m] = V*math.cos(theta) z[m] = V*math.cos(theta)
# imag node voltage
z[m+1] = V*math.sin(theta) z[m+1] = V*math.sin(theta)
H[m][0] = 1; H[m][0] = 1;
H[m][2:] = np.concatenate((-znod.R[i],znod.X[i]), axis=0) H[m][2:] = np.concatenate((-znod.R[i],znod.X[i]),axis=0)
H[m+1][1] = 1 H[m+1][1] = 1
H[m+1][2:] = np.concatenate((-znod.X[i],-znod.R[i]), axis=0) H[m+1][2:] = np.concatenate((-znod.X[i],-znod.R[i]),axis=0)
elif node.type[i] == 'PQ': elif nodes.type[i] == 'PQ':
to = np.where(branch.end==k) to = np.where(branches.end==k)
fr = np.where(branch.start==k) fr = np.where(branches.start==k)
to1 = to[0]+2 to1 = to[0]+2
to2 = to1+branch.num to2 = to1+branches.num
fr1 = fr[0]+2 fr1 = fr[0]+2
fr2 = fr1+branch.num fr2 = fr1+branches.num
H[m][to1] = 1 H[m][to1] = 1
H[m+1][to2] = 1 H[m+1][to2] = 1
H[m][fr1] = -1 H[m][fr1] = -1
H[m+1][fr2] = -1 H[m+1][fr2] = -1
elif node.type[i] == 'PV': elif nodes.type[i] == 'PV':
z[m+1] = node.pwr_flow_values[2][i] z[m+1] = nodes.pwr_flow_values[2][i]
to = np.where(branch.end==k) to = np.where(branches.end==k)
fr = np.where(branch.start==k) fr = np.where(branches.start==k)
to1 = to[0]+2 to1 = to[0]+2
fr1 = fr[0]+2 fr1 = fr[0]+2
H[m][to1] = 1 H[m][to1] = 1
H[m][fr1] = -1 H[m][fr1] = -1
print(z) diff = 5
print(h) epsilon = 10**(-10)
print(H) V = np.ones(nodes.num) + 1j* np.zeros(nodes.num)
epsilon = 5
V = np.ones(node.num) + 1j* np.zeros(node.num)
num_iter = 0 num_iter = 0
State = np.zeros(2*branch.num) State = np.zeros(2*branches.num)
State = np.concatenate( (np.array([1,0]),State), axis=0 ) State = np.concatenate((np.array([1,0]),State),axis=0)
while epsilon > 10**(-10): while diff > epsilon:
for k in range(1,node.num+1): for k in range(1,nodes.num+1):
i = k-1 i = k-1
m = 2*i m = 2*i
t = node.type[i] t = nodes.type[i]
if t == 'slack': if t == 'slack':
h[m] = np.inner(H[m],State) h[m] = np.inner(H[m],State)
h[m+1] = np.inner(H[m+1],State) h[m+1] = np.inner(H[m+1],State)
elif t == 'PQ': elif t == 'PQ':
z[m] = (node.pwr_flow_values[1][i]*V[i].real + node.pwr_flow_values[2][i]*V[i].imag)/(np.abs(V[i])**2) z[m] = (nodes.pwr_flow_values[1][i]*V[i].real + nodes.pwr_flow_values[2][i]*V[i].imag)/(np.abs(V[i])**2)
z[m+1] = (node.pwr_flow_values[1][i]*V[i].imag - node.pwr_flow_values[2][i]*V[i].real)/(np.abs(V[i])**2) z[m+1] = (nodes.pwr_flow_values[1][i]*V[i].imag - nodes.pwr_flow_values[2][i]*V[i].real)/(np.abs(V[i])**2)
h[m] = np.inner(H[m],State) h[m] = np.inner(H[m],State)
h[m+1] = np.inner(H[m+1],State) h[m+1] = np.inner(H[m+1],State)
elif t == 'PV': elif t == 'PV':
z[m] = (node.pwr_flow_values[1][i]*V[i].real + node.pwr_flow_values[2][i]*V[i].imag)/(np.abs(V[i])**2) z[m] = (nodes.pwr_flow_values[1][i]*V[i].real + nodes.pwr_flow_values[2][i]*V[i].imag)/(np.abs(V[i])**2)
h[m] = np.inner(H[m],State) h[m] = np.inner(H[m],State)
h[m+1] = np.abs(V[i]) h[m+1] = np.abs(V[i])
H[m+1][0] = np.cos(np.angle(V[i])) H[m+1][0] = np.cos(np.angle(V[i]))
H[m+1][1] = np.sin(np.angle(V[i])) H[m+1][1] = np.sin(np.angle(V[i]))
idx = np.where(znod.Z[i] != 0) idx = np.where(znod.Z[i] != 0)
idx1 = idx[0]+2 idx1 = idx[0]+2
idx2 = idx1[0]+branch.num idx2 = idx1[0]+branches.num
H[m+1][idx1] = - znod.R[i][idx]*np.cos(np.angle(V[i])) - znod.X[i][idx]*np.sin(np.angle(V[i])) H[m+1][idx1] = - znod.R[i][idx]*np.cos(np.angle(V[i])) - znod.X[i][idx]*np.sin(np.angle(V[i]))
H[m+1][idx2] = - znod.R[i][idx]*np.sin(np.angle(V[i])) - znod.X[i][idx]*np.cos(np.angle(V[i])) H[m+1][idx2] = - znod.R[i][idx]*np.sin(np.angle(V[i])) - znod.X[i][idx]*np.cos(np.angle(V[i]))
r = np.subtract(z,h) r = np.subtract(z,h)
Hinv = np.linalg.inv(H) Hinv = np.linalg.inv(H)
Delta_State = np.inner(Hinv,r) Delta_State = np.inner(Hinv,r)
State = State + Delta_State State = State + Delta_State
epsilon = np.amax(np.absolute(Delta_State)) diff = np.amax(np.absolute(Delta_State))
V[0] = State[0] + 1j * State[1] V[0] = State[0] + 1j * State[1]
Ir = State[2:branch.num+2] Ir = State[2:branches.num+2]
Ix = State[branch.num+2:] Ix = State[branches.num+2:]
Irx = Ir + 1j*Ix Irx = Ir + 1j*Ix
Vcomplex = np.inner((V[0].real + 1j*V[0].imag), np.ones(node.num)) Vcomplex = np.inner((V[0].real + 1j*V[0].imag), np.ones(nodes.num))
DeltaV = np.inner(znod.Z,Irx) DeltaV = np.inner(znod.Z,Irx)
V = Vcomplex[0] - DeltaV V = Vcomplex[0] - DeltaV
num_iter = num_iter+1 num_iter = num_iter+1
I = Ir + 1j*Ix I = Ir + 1j*Ix
Iinj_r = np.zeros(node.num) Iinj_r = np.zeros(nodes.num)
Iinj_x = np.zeros(node.num) Iinj_x = np.zeros(nodes.num)
for k in range(1,nodes.num+1):
for k in range(1,node.num+1): to = np.where(branches.end==k)
to = np.where(branch.end==k) fr = np.where(branches.start==k)
fr = np.where(branch.start==k)
Iinj_r[k-1] = np.sum(I[to[0]].real) - np.sum(I[fr[0]].real) Iinj_r[k-1] = np.sum(I[to[0]].real) - np.sum(I[fr[0]].real)
Iinj_x[k-1] = np.sum(I[to[0]].imag) - np.sum(I[fr[0]].imag) Iinj_x[k-1] = np.sum(I[to[0]].imag) - np.sum(I[fr[0]].imag)
...@@ -138,8 +126,8 @@ def BC_power_flow(branch, node): ...@@ -138,8 +126,8 @@ def BC_power_flow(branch, node):
Sinj_rx = np.multiply(V, np.conj(Iinj)) Sinj_rx = np.multiply(V, np.conj(Iinj))
Sinj = np.real(Sinj_rx) + 1j * np.imag(Sinj_rx) Sinj = np.real(Sinj_rx) + 1j * np.imag(Sinj_rx)
S1_rx = np.multiply(V[branch.start-1], np.conj(I)) S1_rx = np.multiply(V[branchs.start-1], np.conj(I))
S2_rx = np.multiply(V[branch.end-1], np.conj(I)) S2_rx = np.multiply(V[branchs.end-1], np.conj(I))
S1 = np.real(S1_rx) + 1j * np.imag(S1_rx) S1 = np.real(S1_rx) + 1j * np.imag(S1_rx)
S2 = np.real(S2_rx) + 1j * np.imag(S2_rx) S2 = np.real(S2_rx) + 1j * np.imag(S2_rx)
......
acs-state-estimation/nv_powerflow.py
View file @ 3db666e1
...@@ -23,140 +23,6 @@ def Ymatrix_calc(branch, node): ...@@ -23,140 +23,6 @@ def Ymatrix_calc(branch, node):
Adjacencies[to].append(fr+1) Adjacencies[to].append(fr+1)
return Ymatrix, Adjacencies return Ymatrix, Adjacencies
def BC_power_flow(branch, node):
"""It performs Power Flow by using rectangular branch current state variables."""
class Znod:
def __init__(self,branch,node):
self.Z = numpy.zeros((node.num-1,branch.num),dtype=numpy.complex)
for k in range(1,node.num):
i = k-1
aa = branch.start[i]-2
if aa == -1:
self.Z[i] = self.Z[i]
else:
self.Z[i] = self.Z[aa]
self.Z[i][i] = branch.Z[i]
zzz = numpy.zeros((1,branch.num))
self.Z = numpy.concatenate((zzz,self.Z), axis=0)
self.R = self.Z.real
self.X = self.Z.imag
znod = Znod(branch,node)
z = numpy.zeros(2*(branch.num+1))
h = numpy.zeros(2*(branch.num+1))
H = numpy.zeros((2*(branch.num+1),2*(branch.num+1)))
for k in range(1,node.num+1):
i = k-1
m = 2*i
t = node.type[i]
if t == 'slack':
V = node.pwr_flow_values[1][i]
theta = node.pwr_flow_values[2][i]
z[m] = V*math.cos(theta)
z[m+1] = V*math.sin(theta)
H[m][0] = 1;
H[m][2:] = numpy.concatenate((-znod.R[i],znod.X[i]),axis=0)
H[m+1][1] = 1
H[m+1][2:] = numpy.concatenate((-znod.X[i],-znod.R[i]),axis=0)
elif t == 'PQ':
to = numpy.where(branch.end==k)
fr = numpy.where(branch.start==k)
to1 = to[0]+2
to2 = to1+branch.num
fr1 = fr[0]+2
fr2 = fr1+branch.num
H[m][to1] = 1
H[m+1][to2] = 1
H[m][fr1] = -1
H[m+1][fr2] = -1
elif t == 'PV':
z[m+1] = node.pwr_flow_values[2][i]
to = numpy.where(branch.end==k)
fr = numpy.where(branch.start==k)
to1 = to[0]+2
fr1 = fr[0]+2
H[m][to1] = 1
H[m][fr1] = -1
epsilon = 5
Vr = numpy.ones(node.num)
Vx = numpy.zeros(node.num)
V = Real_to_all(Vr, Vx)
num_iter = 0
State = numpy.zeros(2*branch.num)
State = numpy.concatenate((numpy.array([1,0]),State),axis=0)
while epsilon>10**(-10):
for k in range(1,node.num+1):
i = k-1
m = 2*i
t = node.type[i]
if t == 'slack':
h[m] = numpy.inner(H[m],State)
h[m+1] = numpy.inner(H[m+1],State)
elif t == 'PQ':
z[m] = (node.pwr_flow_values[1][i]*V.real[i] + node.pwr_flow_values[2][i]*V.imag[i])/(V.mag[i]**2)
z[m+1] = (node.pwr_flow_values[1][i]*V.imag[i] - node.pwr_flow_values[2][i]*V.real[i])/(V.mag[i]**2)
h[m] = numpy.inner(H[m],State)
h[m+1] = numpy.inner(H[m+1],State)
elif t == 'PV':
z[m] = (node.pwr_flow_values[1][i]*V.real[i] + node.pwr_flow_values[2][i]*V.imag[i])/(V.mag[i]**2)
h[m] = numpy.inner(H[m],State)
h[m+1] = V.mag[i]
H[m+1][0] = numpy.cos(V.phase[i])
H[m+1][1] = numpy.sin(V.phase[i])
idx = numpy.where(znod.Z[i] != 0)
idx1 = idx[0]+2
idx2 = idx1[0]+branch.num
H[m+1][idx1] = - znod.R[i][idx]*numpy.cos(V.phase[i]) - znod.X[i][idx]*numpy.sin(V.phase[i])
H[m+1][idx2] = - znod.R[i][idx]*numpy.sin(V.phase[i]) - znod.X[i][idx]*numpy.cos(V.phase[i])
r = numpy.subtract(z,h)
Hinv = numpy.linalg.inv(H)
Delta_State = numpy.inner(Hinv,r)
State = State + Delta_State
epsilon = numpy.amax(numpy.absolute(Delta_State))
V.real[0] = State[0]
V.imag[0] = State[1]
Ir = State[2:branch.num+2]
Ix = State[branch.num+2:]
Irx = Ir + 1j*Ix
Vcomplex = numpy.inner((V.real[0] + 1j*V.imag[0]),numpy.ones(node.num))
DeltaV = numpy.inner(znod.Z,Irx)
V.complex = Vcomplex[0] - DeltaV
V.real = numpy.real(V.complex)
V.imag = numpy.imag(V.complex)
V = Real_to_all(V.real, V.imag)
num_iter = num_iter+1
I = Real_to_all(Ir,Ix)
Iinj_r = numpy.zeros(node.num)
Iinj_x = numpy.zeros(node.num)
for k in range(1,node.num+1):
to = numpy.where(branch.end==k)
fr = numpy.where(branch.start==k)
Iinj_r[k-1] = numpy.sum(I.real[to[0]]) - numpy.sum(I.real[fr[0]])
Iinj_x[k-1] = numpy.sum(I.imag[to[0]]) - numpy.sum(I.imag[fr[0]])
Iinj = Real_to_all(Iinj_r,Iinj_x)
Sinj_rx = numpy.multiply(V.complex,numpy.conj(Iinj.complex))
Sinj = Real_to_all(numpy.real(Sinj_rx),numpy.imag(Sinj_rx))
S1_rx = numpy.multiply(V.complex[branch.start-1],numpy.conj(I.complex))
S2_rx = - numpy.multiply(V.complex[branch.end-1],numpy.conj(I.complex))
S1 = Real_to_all(numpy.real(S1_rx),numpy.imag(S1_rx))
S2 = Real_to_all(numpy.real(S2_rx),numpy.imag(S2_rx))
return V, I, Iinj, S1, S2, Sinj, num_iter
def NV_power_flow(branch, node): def NV_power_flow(branch, node):
"""It performs Power Flow by using rectangular node voltage state variables.""" """It performs Power Flow by using rectangular node voltage state variables."""
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment