Wind farm

In [207]:
import numpy as np
from pydgrid import grid
from pydgrid.plot_bokeh import plot_results

System

In [210]:
p_gen = 2000.0  # kW
q_gen = 0       # kvar

p_statcom = 0.0  # kW
q_statcom = 0.0     # kvar

data = {
    "lines":[
        {"bus_j": "W1mv",  "bus_k": "W2mv",   "code": "mv_al_300", "m":500},
        {"bus_j": "W2mv",  "bus_k": "W3mv",   "code": "mv_al_300", "m":500},
        {"bus_j": "W3mv",  "bus_k": "POImv",  "code": "mv_al_300", "m":500},
        {"bus_j": "POI",  "bus_k": "GRID",  "code": "hv_line", "m":50.0e3},
            ],
    "buses":[
            {"bus": "W1lv",  "pos_x": -1500.0, "pos_y":  200.0, "units": "m", "U_kV":0.69},
            {"bus": "W2lv",  "pos_x": -1000.0, "pos_y":  200.0, "units": "m", "U_kV":0.69},
            {"bus": "W3lv",  "pos_x":  -500.0, "pos_y":  200.0, "units": "m", "U_kV":0.69},
            {"bus": "W1mv",  "pos_x": -1500.0, "pos_y":  180.0, "units": "m", "U_kV":20.0},
            {"bus": "W2mv",  "pos_x": -1000.0, "pos_y":  180.0, "units": "m", "U_kV":20.0},
            {"bus": "W3mv",  "pos_x":  -500.0, "pos_y":  180.0, "units": "m", "U_kV":20.0},
            {"bus": "POImv", "pos_x":     0.0, "pos_y":    0.0, "units": "m", "U_kV":20.0},
            {"bus": "POI",   "pos_x":   100.0, "pos_y":    0.0, "units": "m", "U_kV":66.0},
            {"bus": "GRID",  "pos_x":   500.0, "pos_y":    0.0, "units": "m", "U_kV":66.0},
                    ],
    "transformers":[
                    {"bus_j": "POImv",  "bus_k": "POI",  "S_n_kVA": 10000.0, "U_1_kV":20.0, "U_2_kV":66.0,
                    "R_cc_pu": 0.01, "X_cc_pu":0.08, "connection": "Dyg11_3w",   "conductors_1": 3, "conductors_2": 3},
                    {"bus_j": "W1mv",  "bus_k": "W1lv",  "S_n_kVA": 2500.0, "U_1_kV":20, "U_2_kV":0.69,
                    "R_cc_pu": 0.01, "X_cc_pu":0.06, "connection": "Dyg11_3w",   "conductors_1": 3, "conductors_2": 3},
                    {"bus_j": "W2mv",  "bus_k": "W2lv",  "S_n_kVA": 2500.0, "U_1_kV":20, "U_2_kV":0.69,
                    "R_cc_pu": 0.01, "X_cc_pu":0.06, "connection": "Dyg11_3w",   "conductors_1": 3, "conductors_2": 3},
                    {"bus_j": "W3mv",  "bus_k": "W3lv",  "S_n_kVA": 2500.0, "U_1_kV":20, "U_2_kV":0.69,
                    "R_cc_pu": 0.01, "X_cc_pu":0.06, "connection": "Dyg11_3w",   "conductors_1": 3, "conductors_2": 3},
                    ],
    "grid_formers":[
                    {"bus": "GRID","bus_nodes": [1, 2, 3],
                            "kV": [38.105, 38.105, 38.105], "deg": [30, 150, 270.0]
                    }
                    ],
    "grid_feeders":[{"bus": "W1lv","bus_nodes": [1, 2, 3],"kW": p_gen, "kvar": q_gen},
                    {"bus": "W2lv","bus_nodes": [1, 2, 3],"kW": p_gen, "kvar": q_gen},
                    {"bus": "W3lv","bus_nodes": [1, 2, 3],"kW": p_gen, "kvar": q_gen},
                    {"bus": "POImv","bus_nodes": [1, 2, 3],"kW": p_statcom, "kvar": q_statcom} # STATCOM
    ],
    "groundings":[
        {"bus": "POImv" , "R_gnd":32.0, "X_gnd":0.0, "conductors": 3}
       ],
    "line_codes":
        {
        "mv_al_150":   {"R1":0.262,  "X1":0.118, "C_1_muF":0.250 },
        "mv_al_185":   {"R1":0.209,  "X1":0.113, "C_1_muF":0.281 },
        "mv_al_240":   {"R1":0.161,  "X1":0.109, "C_1_muF":0.301 },
        "mv_al_300":   {"R1":0.128,  "X1":0.105, "C_1_muF":0.340 },
         "hv_line":  {"R1":0.219, "X1":0.365, "R0":0.219,   "X0":0.365}
        }
    }
In [211]:
grid_1 = grid()
grid_1.read(data)  # Load data
grid_1.pf()  # solve power flow
p=plot_results(grid_1)
In [214]:
mon = grid_1.monitor(bus_from='POI',bus_to='GRID')
mon.P

Out[214]:
5910895.785662318
In [215]:
mon.Q

Out[215]:
-490329.45241958584

Short circuits

Three phase at MV side POI

In [4]:
data['shunts'] = [ # three phase fault to ground
                 {"bus": "POImv" , "R":1.0e-8, "X": 0.0, "bus_nodes": [1,2]},
                 {"bus": "POImv" , "R":1.0e-8, "X": 0.0, "bus_nodes": [2,3]},
                 {"bus": "POImv" , "R":1.0e-8, "X": 0.0, "bus_nodes": [3,0]},
                 ]

# powers to zero:
data['grid_feeders'] = [{"bus": "W1lv","bus_nodes": [1, 2, 3],"kW": 0, "kvar": 0},
                        {"bus": "W2lv","bus_nodes": [1, 2, 3],"kW": 0, "kvar": 0},
                        {"bus": "W3lv","bus_nodes": [1, 2, 3],"kW": 0, "kvar": 0},
                        {"bus": "POImv","bus_nodes": [1, 2, 3],"kW":0, "kvar": 0}] # STATCOM

grid_1 = grid()
grid_1.read(data)
grid_1.pf()

p=plot_results(grid_1)
In [ ]:
I_cc = grid_1.transformers[0]['i_1a_m']
print('Three phase short circuit current at POImv = {:0.2f} kA'.format(I_cc/1000))
Three phase short circuit current at POImv = 2.28 kA

Phase-ground W1mv bus

In [ ]:
data['shunts'] = [
                 {"bus": "W1mv" , "R": 1.0e-8, "X": 0.0, "bus_nodes": [1,0]},
                 ]
grid_1 = grid()
grid_1.read(data)  # Load data
grid_1.pf()

p=plot_results(grid_1)
In [ ]:
I_cc = grid_1.transformers[0]['i_1a_m']
print('Phase-ground short circuit current at W1mv = {:0.2f} kA'.format(I_cc/1000))
Phase-ground short circuit current at W1mv = 0.62 kA

Get POI voltage with generators reactive powers

In [ ]:
from scipy import optimize as sopt

data['shunts'] = []

V_ref = 1.0
p_gen = 2.0e3
q_gen = 0

data['grid_feeders'][0]['kW'] = p_gen
data['grid_feeders'][1]['kW'] = p_gen
data['grid_feeders'][2]['kW'] = p_gen

grid_1 = grid()
grid_1.read(data)  # Load data
grid_1.pf()

def residual(x):

    q_gen = x
    data['grid_feeders'][0]['kvar'] = q_gen
    data['grid_feeders'][1]['kvar'] = q_gen
    data['grid_feeders'][2]['kvar'] = q_gen
    grid_1.read(data)  # Load data
    grid_1.pf()  # solve power flow

    V = abs(grid_1.res['POI'].v_ag)/66.0e3*np.sqrt(3)

    return V_ref - V


res = sopt.bisect(residual,-3000.0,3000.0)
res

p=plot_results(grid_1)

Optimization with reactive powers (without STATCOM)

In [ ]:
V_ref = 1.0
p_gen = 3.0e3
q_gen = 0

data['grid_feeders'][0]['kW'] = p_gen
data['grid_feeders'][1]['kW'] = p_gen
data['grid_feeders'][2]['kW'] = p_gen

grid_1 = grid()
grid_1.read(data)  # Load data
grid_1.pf()


def obj(x):

    data['grid_feeders'][0]['kvar'] = x[0]
    data['grid_feeders'][1]['kvar'] = x[1]
    data['grid_feeders'][2]['kvar'] = x[2]
    grid_1.read(data)  # Load data
    grid_1.pf()  # solve power flow
    mon = grid_1.monitor(bus_from='POI', bus_to='GRID')

    P_loss = p_gen*3*1000 - mon.P

    return P_loss


res = sopt.minimize(obj,[0,0,0], method='SLSQP')
print(res)

p=plot_results(grid_1)
     fun: 192203.64692081884
     jac: array([0.875, 0.875, 1.375])
 message: 'Optimization terminated successfully.'
    nfev: 61
     nit: 8
    njev: 8
  status: 0
 success: True
       x: array([496.46556547, 496.46556547, 466.3287252 ])
In [ ]:
mon = grid_1.monitor(bus_from='POI', bus_to='GRID')
print('POI active power = {:0.3f} MW'.format(mon.P/1e6))
print('POI reactive power = {:0.2f} Mvar'.format(mon.Q/1e6))
POI active power = 8.808 MW
POI reactive power = 0.32 Mvar

Optimization with reactive powers (with STATCOM)

In [ ]:
V_ref = 1.0
p_gen = 3.0e3
q_gen = 0

data['grid_feeders'][0]['kW'] = p_gen
data['grid_feeders'][1]['kW'] = p_gen
data['grid_feeders'][2]['kW'] = p_gen

grid_1 = grid()
grid_1.read(data)  # Load data
grid_1.pf()


def obj(x):


    data['grid_feeders'][0]['kvar'] = x[0]
    data['grid_feeders'][1]['kvar'] = x[1]
    data['grid_feeders'][2]['kvar'] = x[2]
    data['grid_feeders'][3]['kvar'] = x[3]
    grid_1.read(data)  # Load data
    grid_1.pf()  # solve power flow
    mon = grid_1.monitor(bus_from='POI', bus_to='GRID')

    P_loss = p_gen*3*1000 - mon.P

    return P_loss

def const1(x):
    return 1.02-abs(grid_1.monitor(bus_from='POI', bus_to='GRID').V_a)/38.105/1000

def const2(x):
    return -(0.98-abs(grid_1.monitor(bus_from='POI', bus_to='GRID').V_a)/38.105/1000)


res = sopt.minimize(obj,[0,0,0,0], method='COBYLA',
                    constraints=[{'type':'ineq','fun':const1},{'type':'ineq','fun':const2}]
                   )
print(res)

p=plot_results(grid_1)
     fun: 192924.1092131082
   maxcv: 0.0
 message: 'Optimization terminated successfully.'
    nfev: 644
  status: 1
 success: True
       x: array([201.34719703, 200.97303112, 200.48882721, 286.99544409])
In [ ]:
mon = grid_1.monitor(bus_from='POI', bus_to='GRID')
print('POI active power = {:0.3f} MW'.format(mon.P/1e6))
print('POI reactive power = {:0.2f} Mvar'.format(mon.Q/1e6))
POI active power = 8.807 MW
POI reactive power = -0.25 Mvar

Run this notebook in binder:

Documentation Status