2824 lines
126 KiB
Plaintext
2824 lines
126 KiB
Plaintext
|
# HERE WE ARE CONSIDERING DIFFEREN TIME of arival and depature of EV,
|
||
|
|
||
|
import itertools
|
||
|
import pyomo.environ as py
|
||
|
import pandas as pd
|
||
|
import numpy as np
|
||
|
from matplotlib import pyplot as plt
|
||
|
from matplotlib import rc
|
||
|
|
||
|
# to ignore the warning for dividing by 0
|
||
|
np.seterr(divide='ignore', invalid='ignore')
|
||
|
m1 = py.ConcreteModel()
|
||
|
|
||
|
# <editor-fold desc="ReadingTHE OPF DATA FROM EXCEL for HOME">
|
||
|
sheet_data = pd.read_excel("input_data.xlsx", sheet_name='data')
|
||
|
sheet_tr = pd.read_excel("input_data.xlsx", sheet_name='Ambient_Temp')
|
||
|
sheet_EV = pd.read_excel("input_data.xlsx", sheet_name='EV_100')
|
||
|
# sheet_EV = pd.read_excel("input_data.xlsx", sheet_name='EV_10')
|
||
|
# sheet_EV = pd.read_excel("input_data.xlsx", sheet_name='EV_50')
|
||
|
# sheet_EV = pd.read_excel("input_data.xlsx", sheet_name='EV_80')
|
||
|
sheet_EWH = pd.read_excel("input_data.xlsx", sheet_name='EWH')
|
||
|
sheet_wateruse = pd.read_excel("input_data.xlsx", sheet_name='Water_Use')
|
||
|
sheet_home = pd.read_excel("input_data.xlsx", sheet_name='Home')
|
||
|
sheet_ESS = pd.read_excel("input_data.xlsx", sheet_name='ESS_100')
|
||
|
# sheet_ESS = pd.read_excel("input_data.xlsx", sheet_name='ESS_10')
|
||
|
# sheet_ESS = pd.read_excel("input_data.xlsx", sheet_name='ESS_50')
|
||
|
# sheet_ESS = pd.read_excel("input_data.xlsx", sheet_name='ESS_80')
|
||
|
sheet_PV = pd.read_excel("input_data.xlsx", sheet_name='PV_100')
|
||
|
# sheet_PV = pd.read_excel("input_data.xlsx", sheet_name='PV_10')
|
||
|
# sheet_PV = pd.read_excel("input_data.xlsx", sheet_name='PV_50')
|
||
|
# sheet_PV = pd.read_excel("input_data.xlsx", sheet_name='PV_80')
|
||
|
sheet_load = pd.read_excel("input_data.xlsx", sheet_name='Base_Load')
|
||
|
sheet_bus = pd.read_excel("input_data.xlsx", sheet_name='Bus')
|
||
|
sheet_bus1 = pd.read_excel("input_data.xlsx", sheet_name='bus1')
|
||
|
# sheet_bus1 = pd.read_excel("input_data.xlsx", sheet_name='bus1')
|
||
|
sheet_sumdict = pd.read_excel("input_data.xlsx", sheet_name='sum_dict')
|
||
|
sheet_ambient = pd.read_excel("input_data.xlsx", sheet_name='Ambient_Temp')
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
|
||
|
# <editor-fold desc="Home Model for Minimizing Cost">
|
||
|
|
||
|
# <editor-fold desc="Set and Parameters of HOME">
|
||
|
m1.t = py.Set(initialize=sheet_data.time, ordered=True, doc='time period')
|
||
|
m1.b = py.Set(initialize=sheet_bus.bus, doc='bus')
|
||
|
m1.b1 = py.Set(initialize=sheet_bus1.bus1, doc='bus')
|
||
|
# 6:30 = 39
|
||
|
m1.h = py.Set(initialize=sheet_home.home)
|
||
|
|
||
|
m1.cost = py.Param(m1.h, mutable=True)
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="base load and other Parameters of Home">
|
||
|
# Base load : uncontrollable loads
|
||
|
m1.base_load = py.Param(m1.h, m1.t,
|
||
|
initialize=dict(zip(list(itertools.product(m1.h.data(), m1.t.data())),
|
||
|
[i for x in sheet_load.columns if "load" in x for i in
|
||
|
sheet_load.loc[:, x].values])), doc="Base Load")
|
||
|
|
||
|
m1.PV = py.Param(m1.h, m1.t, initialize=dict(zip(list(itertools.product(m1.h.data(), m1.t.data())),
|
||
|
[i for x in sheet_PV.columns if "PV" in x for i in
|
||
|
sheet_PV.loc[:, x].values])), doc="PV Production")
|
||
|
m1.max = py.Param(initialize=1000000, doc='maximum value for selling and buying power')
|
||
|
# m1.Transformer_Limit = py.Param(initialize=75, doc='Power limit of transformer')
|
||
|
|
||
|
m1.UDC = py.Param(m1.t, mutable=True, initialize=80, doc='Uniform Distribution Capacity of transformer')
|
||
|
m1.Buy_Price_p2p = py.Param(m1.h, m1.t, initialize=0, mutable=True)
|
||
|
m1.combine = py.Param(m1.t, initialize=0, mutable=True)
|
||
|
m1.Sell_Price_p2p = py.Param(m1.h, m1.t, initialize=0, mutable=True)
|
||
|
m1.P_Buy_extra = py.Param(m1.h, m1.t, initialize=0)
|
||
|
m1.P_Buy_extra_1 = py.Param(m1.h, m1.t, initialize=0)
|
||
|
m1.P_sell_extra = py.Param(m1.h, m1.t, initialize=0)
|
||
|
m1.P_sell_extra_1 = py.Param(m1.h, m1.t, initialize=0)
|
||
|
m1.Sell_price_p2p = py.Param(m1.h, m1.t, initialize=0)
|
||
|
# Calculating the uniform distribution capacity
|
||
|
# m1.Transformer_Limit = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Transformer_limit)),
|
||
|
# doc='Maximum power limit on transformer')
|
||
|
|
||
|
# for i in m1.t:
|
||
|
# m1.UDC[i] = m1.Transformer_Limit[i] / m1.h[-1]
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Electric Vehicle Parameters">
|
||
|
# m1.alpha = py.Param(initialize=37, doc='time of arrival ') # (40 = 7:00 pm) 37 = 6 pm
|
||
|
m1.alpha = py.Param(m1.h, initialize=dict(zip(list(itertools.product(m1.h.data())),
|
||
|
[i for x in sheet_EV.columns if "arrival" in x for i
|
||
|
in
|
||
|
sheet_EV.loc[:, x].values])), doc="Arrival time of EV")
|
||
|
|
||
|
m1.beta = py.Param(m1.h, initialize=dict(zip(list(itertools.product(m1.h.data())),
|
||
|
[i for x in sheet_EV.columns if "departure" in x for i
|
||
|
in
|
||
|
sheet_EV.loc[:, x].values])), doc="departure time of EV")
|
||
|
# m1.beta = py.Param(initialize=92, doc='time of departure ') # (92 = 8:00 am)
|
||
|
m1.t_final = py.Param(initialize=96, doc='final time of the day ') # 96
|
||
|
m1.second = py.Param(initialize=2, doc='second time of the day ')
|
||
|
m1.EV_Dur = py.Set(m1.h, initialize=lambda m, h: range(m.alpha[h], m.beta[h] + 1), doc='time interval for EV')
|
||
|
# m1.EV_Dur = py.Set(m1.h, initialize=py.RangeSet((m1.alpha[h].value, m1.beta[h].value) for h in m1.h),
|
||
|
# doc='time interval for EV')
|
||
|
# m1.EV_Dur = py.Set(initialize=((h, t) for h in m1.h for t in range(m1.alpha[h], m1.beta[h]+1)), doc='time interval for EV')
|
||
|
|
||
|
# m1.EV_Dur = py.Set(initialize={(i,j) for i in m1.h for j in range(m1.alpha[i].value, m1.beta[i].value + 1)},
|
||
|
# dimen=2,
|
||
|
# doc='time interval for EV')
|
||
|
|
||
|
|
||
|
m1.EV_Dur1 = py.Set(m1.h, initialize=lambda m, h: range(m.alpha[h] + 1, m.beta[h] + 1), doc='time interval for EV')
|
||
|
#
|
||
|
#
|
||
|
# def init_EV_Dur1(m, h):
|
||
|
# return range(m.alpha[h] + 1, m.beta[h] + 1)
|
||
|
#
|
||
|
# m1.EV_Dur1 = py.Set(m1.h, initialize=init_EV_Dur1, ordered=True, dimen=1, within=OrderedSimpleSet, doc='time interval for EV')
|
||
|
|
||
|
|
||
|
# m1.EV_Dur1 = py.Set(m1.h,initialize=py.RangeSet(m1.alpha[h].value + 1, m1.beta[h].value),
|
||
|
# doc='time interval for EV after arrival')
|
||
|
m1.t_second = py.Set(initialize=py.RangeSet(m1.second.value, m1.t_final.value),
|
||
|
doc='time greater then 1')
|
||
|
m1.t_first = py.Param(initialize=1, doc='first time of the day ')
|
||
|
# m1.ist_interval = py.Set(initialize=py.RangeSet(m1.t_first.value, m1.alpha - 1),
|
||
|
# doc='time period before the arrival time')
|
||
|
|
||
|
m1.ist_interval = py.Set(m1.h, initialize=lambda m, h: range(m1.t_first.value, m.alpha[h]), doc='time interval for EV')
|
||
|
# m1.second_interval = py.Set(initialize=py.RangeSet(m1.beta.value + 1, m1.t_final.value),
|
||
|
# doc='time interval after departure')
|
||
|
|
||
|
m1.second_interval = py.Set(m1.h, initialize=lambda m, h: range(m.beta[h] + 1, m1.t_final.value + 1),
|
||
|
doc='time interval after departure')
|
||
|
|
||
|
m1.ChargeRate_EV = py.Param(m1.h, initialize=dict(zip(list(itertools.product(m1.h.data())),
|
||
|
[i for x in sheet_EV.columns if "charging_rate" in x for i
|
||
|
in
|
||
|
sheet_EV.loc[:, x].values])), doc="charging rate of EV")
|
||
|
|
||
|
m1.DischRate_EV = py.Param(m1.h, initialize=dict(zip(list(itertools.product(m1.h.data())),
|
||
|
[i for x in sheet_EV.columns if "rate_of_discharging" in x
|
||
|
for i in
|
||
|
sheet_EV.loc[:, x].values])), doc="Discharging rate of EV")
|
||
|
|
||
|
m1.Ch_Effc_EV = py.Param(m1.h, initialize=dict(zip(list(itertools.product(m1.h.data())),
|
||
|
[i for x in sheet_EV.columns if "charging_efficiency" in x for
|
||
|
i in
|
||
|
sheet_EV.loc[:, x].values])), doc="charging efficency of EV")
|
||
|
|
||
|
m1.DischEffc_EV = py.Param(m1.h, initialize=dict(zip(list(itertools.product(m1.h.data())),
|
||
|
[i for x in sheet_EV.columns if
|
||
|
"efficiency_of_discharging" in x for i
|
||
|
in sheet_EV.loc[:, x].values])),
|
||
|
doc="Discharging efficency of EV")
|
||
|
|
||
|
m1.Cap_EV = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EV.columns if "capacity" in x for i in
|
||
|
sheet_EV.loc[:, x].values])), doc="Capacity of EV")
|
||
|
m1.End_percentage_EV = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EV.columns if "end" in x
|
||
|
for i in
|
||
|
sheet_EV.loc[:, x].values])),
|
||
|
doc="Departure energy of EV")
|
||
|
|
||
|
m1.In_Percentage_EV = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EV.columns if "initial" in x
|
||
|
for i in
|
||
|
sheet_EV.loc[:, x].values])),
|
||
|
doc="Initial energy of EV")
|
||
|
|
||
|
m1.Energy_EV_dep = py.Param(m1.h, initialize=dict(
|
||
|
zip(m1.h.data(), np.array(m1.Cap_EV.values()) * m1.End_percentage_EV.values())), # just changing np>py
|
||
|
doc="Departure temperature of EV")
|
||
|
|
||
|
m1.Energy_EV_In = py.Param(m1.h, initialize=dict(
|
||
|
zip(m1.h.data(), np.array(m1.Cap_EV.values()) * m1.In_Percentage_EV.values())),
|
||
|
doc="Initial energy of EV")
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Energy Storage System Parameter">
|
||
|
|
||
|
m1.ChargeRate_ESS = py.Param(m1.h, initialize=dict(zip(list(itertools.product(m1.h.data())),
|
||
|
[i for x in sheet_ESS.columns if "charging_rate" in x for i
|
||
|
in
|
||
|
sheet_ESS.loc[:, x].values])), doc="charging rate of ESS")
|
||
|
|
||
|
m1.DischRate_ESS = py.Param(m1.h, initialize=dict(zip(list(itertools.product(m1.h.data())),
|
||
|
[i for x in sheet_ESS.columns if "rate_of_discharging" in x
|
||
|
for
|
||
|
i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="Discharging rate of ESS")
|
||
|
|
||
|
# m1.ChargeRate_ESS = py.Param(initialize=float(sheet_ESS.charging_rate1), doc='Charging rate of ESS ')
|
||
|
# m1.DischRate_ESS = py.Param(initialize=float(sheet_ESS.discharging_rate1),
|
||
|
# doc='Discharging rate of ESS ')
|
||
|
m1.Cap_ESS = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_ESS.columns if "capacity" in x for i in
|
||
|
sheet_ESS.loc[:, x].values])), doc="Capacity of ESS")
|
||
|
|
||
|
m1.End_percentage_ESS = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_ESS.columns if "end" in x
|
||
|
for i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="Departure energy of ESS")
|
||
|
|
||
|
m1.In_Percentage_ESS = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_ESS.columns if
|
||
|
"initial" in x for i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="Initial energy of ESS")
|
||
|
|
||
|
m1.Ch_Effc_ESS = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_ESS.columns if "charging_efficiency" in x
|
||
|
for i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="charging efficiency of ESS")
|
||
|
|
||
|
m1.DischEffc_ESS = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_ESS.columns if
|
||
|
"efficiency_of_dicharging" in x for
|
||
|
i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="Discharging efficiency of ESS")
|
||
|
|
||
|
m1.Energy_ESS_In = py.Param(m1.h, initialize=dict(
|
||
|
zip(m1.h.data(), np.array(m1.Cap_ESS.values()) * m1.In_Percentage_ESS.values())),
|
||
|
doc="Initial energy of ESS")
|
||
|
|
||
|
m1.End_En_ESS = py.Param(m1.h, initialize=dict(
|
||
|
zip(m1.h.data(), np.array(m1.Cap_ESS.values()) * m1.End_percentage_ESS.values())),
|
||
|
doc="Departure energy of ESS")
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Electric Water Heater Parameter">
|
||
|
|
||
|
m1.tetta_low = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "tetta_low" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])),
|
||
|
doc="Lower bound of the water temperature")
|
||
|
m1.tetta_up = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "tetta_up" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])),
|
||
|
doc="Upper bound of the water temperature")
|
||
|
m1.tetta_amb_int = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "tetta_amb_init" in x for i
|
||
|
in
|
||
|
sheet_EWH.loc[:, x].values])),
|
||
|
doc="Initial ambient temperature")
|
||
|
|
||
|
m1.tetta_amb = py.Param(m1.h, m1.t,
|
||
|
initialize=dict(zip(list(itertools.product(m1.h.data(), m1.t.data())),
|
||
|
[i for x in sheet_ambient.columns if "celsius" in x for i in
|
||
|
sheet_ambient.loc[:, x].values])), doc="Outdoor temperature")
|
||
|
|
||
|
m1.Q = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "capacity" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Power of the EWH")
|
||
|
|
||
|
m1.R = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "thermal_resistance" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Thermal resistance")
|
||
|
m1.C = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "thermal_capacitance" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Thermal resistance")
|
||
|
|
||
|
m1.M = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "water_cap" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Water Capacity (L)")
|
||
|
|
||
|
m1.tetta_EWH_int = py.Param(m1.h, initialize=dict(zip((m1.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "tetta_wat_init" in x for i
|
||
|
in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Initial temperature")
|
||
|
|
||
|
m1.water_use = py.Param(m1.h, m1.t,
|
||
|
initialize=dict(zip(list(itertools.product(m1.h.data(), m1.t.data())),
|
||
|
[i for x in sheet_wateruse.columns if "Litre" in x for i in
|
||
|
sheet_wateruse.loc[:, x].values])), doc="Hot Water usage")
|
||
|
|
||
|
# </editor-fold2
|
||
|
|
||
|
# <editor-fold desc="Price, Base load, and time duration parameters">
|
||
|
# m1.Buy_price = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_RTP)), doc='Buying Price')
|
||
|
# m1.Buy_price = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_TOU)), doc='Buying Price')
|
||
|
# m1.Buy_price = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_TOU_winter)), doc='Buying Price')
|
||
|
# m1.Buy_price = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_flex_D)), doc='Buying Price')
|
||
|
m1.Buy_price = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.TOU_PGE)), doc='Buying Price')
|
||
|
# m1.Buy_price = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.TOU_Texas_RTP)), doc='Buying Price')
|
||
|
# m1.Buy_price_Peer = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_TOU)), doc='Buying Price')
|
||
|
# m1.Buy_price_Peer = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_Peer)),
|
||
|
# doc='Buying Price of peer')
|
||
|
# m1.Sell_price_Peer = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Sell_price_new)),
|
||
|
# doc='Selling Price of peer')
|
||
|
m1.Sell_price = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.FiT)),
|
||
|
doc='Selling Price') # Fit is TX US
|
||
|
# m1.Sell_price_Peer = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.fixed_incentive)), doc='Selling Price')
|
||
|
# Time duration: we took 15 mint granularity so in one hour it will be 1/4
|
||
|
m1.time_d = py.Param(initialize=(1 / 4), doc='time duration ') # 1/4
|
||
|
m1.Export_Percent = py.Param(initialize=1, doc='pecentage export ') # 0.814
|
||
|
|
||
|
m1.cost = py.Param(m1.h, mutable=True)
|
||
|
# </editor-fold>
|
||
|
|
||
|
for h in m1.h:
|
||
|
for t in m1.t:
|
||
|
m1.Buy_Price_p2p[h, t] = m1.Buy_price[t]
|
||
|
for h in m1.h:
|
||
|
for t in m1.t:
|
||
|
m1.Sell_Price_p2p[h, t] = m1.Sell_price[t]
|
||
|
# m1.max_limit = py.Param(m1.t, initialize=75, doc='Maximum power limit on transformer ')
|
||
|
# m1.DPT = py.Param(initialize=0.069, doc='Daily power tariff ') # .069
|
||
|
# </editor-fold>
|
||
|
|
||
|
print('Code Starts for HEMS for cost')
|
||
|
# Variable
|
||
|
|
||
|
# <editor-fold desc="buying and selling variables">
|
||
|
m1.P_Buy_Grid = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
|
||
|
m1.E_Buy_Grid = py.Param(m1.h, mutable=True)
|
||
|
m1.P_Buy_Total = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
# m1.P_Buy_Peer = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.P_Sell_Grid = py.Var(m1.h, m1.t, bounds=(0, 1000000))
|
||
|
m1.P_Sell_Total = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
# m1.P_Sell_Peer = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
|
||
|
m1.S_P_sell = py.Var(m1.h, m1.t, within=py.Binary)
|
||
|
m1.S_P_buy = py.Var(m1.h, m1.t, within=py.Binary)
|
||
|
|
||
|
# m1.pflex = py.Var(m1.b, m1.t, doc='active Flexibility')
|
||
|
# m1.qflex = py.Var(m1.b, m1.t,
|
||
|
# doc='reactive Flexibility') # i comment this because i use m1.pflex variable * 0.48
|
||
|
|
||
|
# Tranformal variable
|
||
|
m1.Power_T_P = py.Var(m1.t, bounds=(0, None), doc="power from transformer to Peer")
|
||
|
m1.Power_Total_power = py.Param(m1.t, mutable=True, doc="power from transformer to Peer")
|
||
|
m1.HST = py.Var(m1.t, bounds=(0, None), doc="HST")
|
||
|
m1.DTOT = py.Var(m1.t, bounds=(0, None), doc="HST")
|
||
|
m1.DHST = py.Var(m1.t, bounds=(0, None), doc="HST")
|
||
|
m1.DTOTU = py.Var(m1.t, bounds=(0, None), doc="HST")
|
||
|
m1.DHSTU = py.Var(m1.t, bounds=(0, None), doc="HST")
|
||
|
m1.Power_P_T = py.Var(m1.t, bounds=(0, None), doc="power from Peer to Transformer")
|
||
|
m1.S_P_T = py.Var(m1.h, m1.t, within=py.Binary)
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Electric vehicle Variables">
|
||
|
m1.P_EV_Charge = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.P_EV_Disch = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.P_EV_Disch_Home = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.P_EV_Disch_Grid = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.S_EV_Charge = py.Var(m1.h, m1.t, within=py.Binary)
|
||
|
m1.S_EV_Disch = py.Var(m1.h, m1.t, within=py.Binary)
|
||
|
m1.Energy_EV = py.Var(m1.h, m1.t, bounds=(0, None)) # SOC of the EV
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Energy storage system Variables">
|
||
|
m1.Energy_ESS = py.Var(m1.h, m1.t, bounds=(0, None)) # SOC of the ESS
|
||
|
m1.S_ESS_Charge = py.Var(m1.h, m1.t, within=py.Binary)
|
||
|
m1.P_ESS_Disch = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.P_ESS_Charge = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.P_ESS_Charge_Grid = py.Var(m1.h, m1.t, bounds=(0, None), doc='ESS charging from the Grid')
|
||
|
m1.P_ESS_Disch_Home = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.P_ESS_Disch_Grid = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.S_ESS_Disch = py.Var(m1.h, m1.t, within=py.Binary)
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="PV variables">
|
||
|
m1.PV_Home = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.PV_Grid = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
m1.PV_Battery = py.Var(m1.h, m1.t, bounds=(0, None))
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="EWH Variables">
|
||
|
m1.tetta_EWH_wat = py.Var(m1.h, m1.t)
|
||
|
m1.S_EWH = py.Var(m1.h, m1.t, within=py.Binary)
|
||
|
m1.P_EWH = py.Var(m1.h, m1.t, doc='power of EWH')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
# Constraints
|
||
|
|
||
|
|
||
|
# <editor-fold desc="Buying and Selling Power">
|
||
|
def Power_buy(m1, h, i):
|
||
|
return m1.P_Buy_Grid[h, i] == m1.base_load[h, i] + m1.Q[h] * m1.S_EWH[h, i] + m1.P_EV_Charge[h, i] - \
|
||
|
m1.P_EV_Disch_Home[h, i] + m1.P_ESS_Charge[h, i] - \
|
||
|
m1.P_ESS_Disch_Home[h, i] - m1.PV_Home[h, i] - m1.PV_Battery[h, i] #
|
||
|
|
||
|
|
||
|
m1.Const_1 = py.Constraint(m1.h, m1.t, rule=Power_buy, doc='Power buy from the Grid')
|
||
|
|
||
|
|
||
|
def Power_buy2(m1, h, i):
|
||
|
return m1.P_Buy_Grid[h, i] <= m1.max * m1.S_P_buy[h, i]
|
||
|
|
||
|
|
||
|
m1.Const_1a = py.Constraint(m1.h, m1.t, rule=Power_buy2,
|
||
|
doc='removing the nonlinearity in the objective fucntion')
|
||
|
|
||
|
|
||
|
def Power_sell1(m1, h, i):
|
||
|
return m1.P_Sell_Grid[h, i] == m1.P_EV_Disch_Grid[h, i] + m1.PV_Grid[h, i] + m1.P_ESS_Disch_Grid[h, i]
|
||
|
|
||
|
|
||
|
m1.Const_2 = py.Constraint(m1.h, m1.t, rule=Power_sell1, doc='Power sell to the Grid')
|
||
|
|
||
|
|
||
|
def Power_sell2(m1, h, i):
|
||
|
# return m1.P_Sell_Grid[h, i] <= m1.max.value * (1 - m1.S_P_buy[h, i])
|
||
|
return m1.P_Sell_Grid[h, i] <= m1.max.value * (1 - m1.S_P_buy[h, i])
|
||
|
|
||
|
|
||
|
m1.Const_2a = py.Constraint(m1.h, m1.t, rule=Power_sell2, doc='Power sell to the Grid')
|
||
|
|
||
|
|
||
|
# put some percentage to export power the grid
|
||
|
def percentage(m1, h, i):
|
||
|
return sum(m1.P_Sell_Grid[h, i] for i in m1.t for h in m1.h) >= m1.Export_Percent.value * sum(
|
||
|
m1.P_EV_Disch_Home[h, i] + m1.P_ESS_Disch_Home[h, i]
|
||
|
+ m1.PV_Battery[h, i] + m1.PV_Home[h, i] for i in m1.t for h in m1.h)
|
||
|
|
||
|
# return sum(m1.P_Sell_Grid[y, i] for i in m1.t for y in m1.h) >= m1.Export_Percent.value * sum(
|
||
|
# m1.P_EV_Disch_Home[y, i] + m1.P_ESS_Disch_Home[y, i]
|
||
|
# + m1.PV_Battery[y, i] + m1.PV_Home[y, i] for i in m1.t for y in m1.h)
|
||
|
|
||
|
|
||
|
m1.connn1 = py.Constraint(m1.h, m1.t, rule=percentage, doc='Export constraint')
|
||
|
|
||
|
|
||
|
#
|
||
|
# def Status_Power(m1, h, i):
|
||
|
# return m1.S_P_buy[h, i] + m1.S_P_sell[h, i] <= 1
|
||
|
#
|
||
|
#
|
||
|
# m1.Const_3 = py.Constraint(m1.h, m1.t, rule=Status_Power,
|
||
|
# doc='Buying and selling power will not occur at same time')
|
||
|
#
|
||
|
|
||
|
# def Power_total_buy(m1, h, i):
|
||
|
# return m1.P_Buy_Total[h, i] == m1.P_Buy_Grid[h, i]
|
||
|
# # return m1.P_Buy_Total[h, i] == m1.P_Buy_Grid[h, i] + m1.P_Buy_Peer[h, i]
|
||
|
#
|
||
|
#
|
||
|
# m1.Const_3a = py.Constraint(m1.h, m1.t, rule=Power_total_buy,
|
||
|
# doc='Total Buying power equal to power buy from grid and form '
|
||
|
# 'other peers')
|
||
|
|
||
|
|
||
|
#
|
||
|
# def Power_total_Sell(m1, h, i):
|
||
|
# return m1.P_Sell_Total[h, i] == m1.P_Sell_Grid[h, i]
|
||
|
# # return m1.P_Sell_Total[h, i] == m1.P_Sell_Grid[h, i] + m1.P_Sell_Peer[h, i]
|
||
|
#
|
||
|
#
|
||
|
# m1.Const_3b = py.Constraint(m1.h, m1.t, rule=Power_total_Sell,
|
||
|
# doc='Total Selling power equal to power sell to grid and to '
|
||
|
# 'other peers')
|
||
|
|
||
|
|
||
|
#
|
||
|
# def Power_total(m1, h, i):
|
||
|
# return sum(m1.P_Buy_Peer[h, i] for h in m1.h) == sum(m1.P_Sell_Peer[h, i] for h in m1.h)
|
||
|
#
|
||
|
#
|
||
|
# m1.Const_3c = py.Constraint(m1.h, m1.t, rule=Power_total,
|
||
|
# doc='Total Buying power equal to total Selling power')
|
||
|
|
||
|
|
||
|
def Power_T_P(m1, h, i):
|
||
|
return m1.Power_T_P[i] == sum(m1.P_Buy_Grid[h, i] for h in m1.h)
|
||
|
|
||
|
|
||
|
#
|
||
|
m1.Const_3d = py.Constraint(m1.h, m1.t, rule=Power_T_P,
|
||
|
doc='Power bought from the grid is equal to the Power of transformer to grid direction')
|
||
|
|
||
|
|
||
|
def Power_P_T(m1, h, i):
|
||
|
return m1.Power_P_T[i] == sum(m1.P_Sell_Grid[h, i] for h in m1.h)
|
||
|
|
||
|
|
||
|
m1.Const_3e = py.Constraint(m1.h, m1.t, rule=Power_P_T,
|
||
|
doc='Power sell to the grid is equal to the Power on the transformer from grid direction')
|
||
|
|
||
|
|
||
|
# Transformer limits
|
||
|
|
||
|
# def Power_limit_T_P(m1, h, i):
|
||
|
# return m1.Power_T_P[i] <= m1.max.value * (m1.S_P_T[h, i])
|
||
|
#
|
||
|
#
|
||
|
# m1.Const_3f = py.Constraint(m1.h, m1.t, rule=Power_limit_T_P, doc='Power limit on the transformer')
|
||
|
#
|
||
|
#
|
||
|
# def Power_limit_P_T(m1, h, i):
|
||
|
# return m1.Power_P_T[i] <= m1.max.value * (1 - m1.S_P_T[h, i])
|
||
|
#
|
||
|
#
|
||
|
# m1.Const_3g = py.Constraint(m1.h, m1.t, rule=Power_limit_P_T, doc='Power limit on the transformer')
|
||
|
|
||
|
|
||
|
# def Power_UDC_Limit(m1, h, i):
|
||
|
# return m1.P_Buy_Grid[h, i] <= m1.UDC[i]
|
||
|
#
|
||
|
#
|
||
|
# m1.Const_3h = py.Constraint(m1.h, m1.t, rule=Power_UDC_Limit, doc='UDC limit')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
#
|
||
|
# <editor-fold desc="Electric vehicle Constraints">
|
||
|
def Power_EV_Charge_limit1(m1, h, i):
|
||
|
return m1.P_EV_Charge[h, i] <= m1.ChargeRate_EV[h] * m1.S_EV_Charge[h, i]
|
||
|
|
||
|
|
||
|
m1.Const_EV_1 = py.Constraint(m1.h, m1.t, rule=Power_EV_Charge_limit1,
|
||
|
doc='Charging power of EV (Upper limit)')
|
||
|
|
||
|
|
||
|
def Power_EV_Charge_limit2(m1, h, i):
|
||
|
return m1.P_EV_Charge[h, i] >= 0
|
||
|
|
||
|
|
||
|
m1.Const_EV_2 = py.Constraint(m1.h, m1.t, rule=Power_EV_Charge_limit2,
|
||
|
doc='Charging power of EV (lower limit)')
|
||
|
|
||
|
|
||
|
def Power_EV_Disch_limit1(m1, h, i):
|
||
|
return m1.P_EV_Disch[h, i] <= (m1.DischRate_EV[h] * m1.S_EV_Disch[h, i])
|
||
|
|
||
|
|
||
|
m1.Const_EV_3 = py.Constraint(m1.h, m1.t, rule=Power_EV_Disch_limit1,
|
||
|
doc='Discharging power of EV (Upper limit)')
|
||
|
|
||
|
|
||
|
def Power_EV_Disch_limit2(m1, h, i):
|
||
|
return m1.P_EV_Disch[h, i] >= 0
|
||
|
|
||
|
|
||
|
m1.Const_EV_4 = py.Constraint(m1.h, m1.t, rule=Power_EV_Disch_limit2,
|
||
|
doc='Discharging power of EV (Lower limit)')
|
||
|
|
||
|
|
||
|
def Status_EV(m1, h, i):
|
||
|
return m1.S_EV_Disch[h, i] + m1.S_EV_Charge[h, i] <= 1
|
||
|
|
||
|
|
||
|
m1.Const_EV_5 = py.Constraint(m1.h, m1.t, rule=Status_EV,
|
||
|
doc='Charging and discharging will not occur at same time')
|
||
|
|
||
|
|
||
|
def Power_EV_Disch(m1, h, i):
|
||
|
return m1.P_EV_Disch[h, i] * m1.DischEffc_EV[h] == (
|
||
|
(m1.P_EV_Disch_Home[h, i] + m1.P_EV_Disch_Grid[h, i]))
|
||
|
|
||
|
|
||
|
m1.Const_EV_6 = py.Constraint(m1.h, m1.t, rule=Power_EV_Disch, doc='Discharging power of EV to '
|
||
|
'Home and Grid')
|
||
|
|
||
|
|
||
|
def SoC_EV1(m1, h, i):
|
||
|
return m1.Energy_EV[h, i] == m1.Energy_EV_In[h]
|
||
|
|
||
|
|
||
|
m1.Const_EV_7a = py.Constraint([(h, py.value(m1.alpha[h])) for h in m1.h], rule=SoC_EV1, doc='SoC of the EV at arrival')
|
||
|
|
||
|
|
||
|
# def SoC_EV1(m1, h):
|
||
|
# alpha_h = m1.alpha[h]
|
||
|
# return m1.Energy_EV[h, alpha_h] == m1.Energy_EV_In[h]
|
||
|
#
|
||
|
# m1.Const_EV_7a = py.Constraint(m1.h, rule=SoC_EV1, doc='SoC of the EV at arrival')
|
||
|
|
||
|
def SoC_EV2(m1, h, i):
|
||
|
return m1.Energy_EV[h, i] == m1.Energy_EV[h, i - 1] + (
|
||
|
m1.P_EV_Charge[h, i] * m1.Ch_Effc_EV[h] * m1.DischEffc_EV[h] - m1.P_EV_Disch[h, i]) * m1.time_d
|
||
|
|
||
|
|
||
|
m1.Const_EV_7b = py.ConstraintList()
|
||
|
for h in m1.h:
|
||
|
for i in m1.EV_Dur1[h].value:
|
||
|
m1.Const_EV_7b.add(SoC_EV2(m1, h, i))
|
||
|
|
||
|
|
||
|
def EV_availability1(m1, h, i):
|
||
|
return m1.Energy_EV[h, i] == 0
|
||
|
|
||
|
|
||
|
# m1.Const_EV_8a = py.Constraint(m1.h, m1.ist_interval[h], rule=EV_availability1,doc='SOC available before arrival time')
|
||
|
# m1.Const_EV_8a = py.Constraint([(h, m1.ist_interval[h].value) for h in m1.h], rule=EV_availability1, doc='SOC available before arrival time')
|
||
|
m1.Const_EV_8a = py.ConstraintList()
|
||
|
for h in m1.h:
|
||
|
for i in m1.ist_interval[h].value:
|
||
|
m1.Const_EV_8a.add(EV_availability1(m1, h, i))
|
||
|
|
||
|
|
||
|
def EV_availability2(m1, h, i):
|
||
|
return m1.Energy_EV[h, i] == 0
|
||
|
|
||
|
|
||
|
# m1.Const_EV_8b = py.Constraint([(h, m1.second_interval[h].value) for h in m1.h], rule=EV_availability2, doc='SOC available after departure time')
|
||
|
# m1.Const_EV_8b = py.Constraint(m1.h, m1.second_interval[h], rule=EV_availability2,doc='SOC available after departure time')
|
||
|
m1.Const_EV_8b = py.ConstraintList()
|
||
|
for h in m1.h:
|
||
|
for i in m1.second_interval[h].value:
|
||
|
m1.Const_EV_8b.add(EV_availability2(m1, h, i))
|
||
|
|
||
|
|
||
|
def EV_status_available1(m1, h, i):
|
||
|
return m1.S_EV_Disch[h, i] + m1.S_EV_Charge[h, i] == 0
|
||
|
|
||
|
|
||
|
# m1.Const_EV_9a = py.Constraint(m1.h, m1.ist_interval[h], rule=EV_status_available1, doc='EV availability before arrival time')
|
||
|
# m1.Const_EV_9a = py.Constraint([(h, m1.ist_interval[h].value) for h in m1.h], rule=EV_status_available1, doc='EV availability before arrival time')
|
||
|
m1.Const_EV_9a = py.ConstraintList()
|
||
|
for h in m1.h:
|
||
|
for i in m1.ist_interval[h].value:
|
||
|
m1.Const_EV_9a.add(EV_status_available1(m1, h, i))
|
||
|
|
||
|
|
||
|
def EV_status_available2(m1, h, i):
|
||
|
return m1.S_EV_Disch[h, i] + m1.S_EV_Charge[h, i] == 0
|
||
|
|
||
|
|
||
|
# m1.Const_EV_9b = py.Constraint(m1.h, m1.second_interval[h], rule=EV_status_available2,doc='EV availability after arrival time')
|
||
|
# m1.Const_EV_9b = py.Constraint([(h, m1.second_interval[h].value) for h in m1.h], rule=EV_status_available2, doc='EV availability after arrival time')
|
||
|
m1.Const_EV_9b = py.ConstraintList()
|
||
|
for h in m1.h:
|
||
|
for i in m1.second_interval[h].value:
|
||
|
m1.Const_EV_9b.add(EV_status_available2(m1, h, i))
|
||
|
|
||
|
|
||
|
def EV_SoC_limit1(m1, h, i):
|
||
|
return m1.Energy_EV[h, i] >= 0.2 * m1.Cap_EV[h]
|
||
|
|
||
|
|
||
|
# m1.Const_EV_10 = py.Constraint(m1.h, m1.EV_Dur[h], rule=EV_SoC_limit1, doc='Minimum SoC of EV')
|
||
|
# m1.Const_EV_10 = py.Constraint([(h, m1.EV_Dur[h].value) for h in m1.h], rule=EV_SoC_limit1, doc='Minimum SoC of EV')
|
||
|
m1.Const_EV_10 = py.ConstraintList()
|
||
|
for h in m1.h:
|
||
|
for i in m1.EV_Dur[h].value:
|
||
|
m1.Const_EV_10.add(EV_SoC_limit1(m1, h, i))
|
||
|
|
||
|
|
||
|
def EV_SoC_limit2(m1, h, i):
|
||
|
return m1.Energy_EV[h, i] <= m1.Cap_EV[h]
|
||
|
|
||
|
|
||
|
# m1.Const_EV_11 = py.Constraint(m1.h, m1.t, rule=EV_SoC_limit2, doc='Maximum SoC of EV')
|
||
|
m1.Const_EV_11 = py.ConstraintList()
|
||
|
for h in m1.h:
|
||
|
for i in m1.EV_Dur[h].value:
|
||
|
m1.Const_EV_11.add(EV_SoC_limit2(m1, h, i))
|
||
|
|
||
|
|
||
|
def EV_final_SoC(m1, h, i):
|
||
|
return m1.Energy_EV[h, i] == m1.Energy_EV_dep[h]
|
||
|
|
||
|
|
||
|
# m1.Const_EV_12 = py.Constraint(m1.h, [m1.beta[h]], rule=EV_final_SoC,doc='Final SoC of EV at departure time')
|
||
|
# m1.Const_EV_12 = py.Constraint([(h,py.value([m1.beta[h]])) for h in m1.h], rule=EV_final_SoC,doc='Final SoC of EV at departure time')
|
||
|
m1.Const_EV_12 = py.Constraint([(h, py.value(m1.beta[h])) for h in m1.h], rule=EV_final_SoC,
|
||
|
doc='Final SoC of EV at departure time')
|
||
|
|
||
|
|
||
|
# m1.Const_EV_12 = py.ConstraintList()
|
||
|
# for h in m1.h:
|
||
|
# for i in py.value(m1.beta[h]):
|
||
|
# m1.Const_EV_12.add(EV_final_SoC(m1, h, i))
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Energy Storage System Constraints">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
||
|
def SoC_ESS1(m1, h, i):
|
||
|
return m1.Energy_ESS[h, i] == m1.Energy_ESS_In[h]
|
||
|
|
||
|
|
||
|
m1.Const_ESS_1a = py.Constraint(m1.h, [m1.t_first.value], rule=SoC_ESS1,
|
||
|
doc='SoC of ESS') # need to check wather we need to put in the brackets or not
|
||
|
|
||
|
|
||
|
def SoC_ESS2(m1, h, i):
|
||
|
return m1.Energy_ESS[h, i] == m1.Energy_ESS[h, i - 1] + (
|
||
|
m1.P_ESS_Charge[h, i] * m1.Ch_Effc_ESS[h] - m1.P_ESS_Disch[h, i]) * m1.time_d
|
||
|
|
||
|
|
||
|
m1.Const_ESS_1b = py.Constraint(m1.h, m1.t_second, rule=SoC_ESS2, doc='SoC of ESS')
|
||
|
|
||
|
|
||
|
def Power_ESS_Charge(m1, h, i):
|
||
|
return m1.P_ESS_Charge[h, i] == m1.P_ESS_Charge_Grid[h, i] + m1.PV_Battery[h, i]
|
||
|
|
||
|
|
||
|
m1.Const_ESS_2 = py.Constraint(m1.h, m1.t, rule=Power_ESS_Charge,
|
||
|
doc='Charging power of ESS ')
|
||
|
|
||
|
|
||
|
def Power_ESS_Charge_limit1(m1, h, i):
|
||
|
return (m1.P_ESS_Charge_Grid[h, i] + m1.PV_Battery[h, i]) <= m1.ChargeRate_ESS[h] * m1.S_ESS_Charge[h, i]
|
||
|
|
||
|
|
||
|
m1.Const_ESS_2a = py.Constraint(m1.h, m1.t, rule=Power_ESS_Charge_limit1,
|
||
|
doc='Charging power of ESS (Upper limit)')
|
||
|
|
||
|
|
||
|
def Power_ESS_Charge_limit2(m1, h, i):
|
||
|
return (m1.P_ESS_Charge_Grid[h, i] + m1.PV_Battery[h, i]) >= 0
|
||
|
# return (m1.P_ESS_Charge_Grid[h, i] + m1.PV_Battery[h, i]) >= 0.2 * (m1.ChargeRate_ESS[h] * m1.S_ESS_Charge[h, i])
|
||
|
|
||
|
|
||
|
m1.Const_ESS_2b = py.Constraint(m1.h, m1.t, rule=Power_ESS_Charge_limit2,
|
||
|
doc='Charging power of ESS (Lower limit)')
|
||
|
|
||
|
|
||
|
def Power_ESS_Disch_limit1(m1, h, i):
|
||
|
return (m1.P_ESS_Disch_Home[h, i] + m1.P_ESS_Disch_Grid[h, i]) <= (m1.DischRate_ESS[h] * m1.S_ESS_Disch[h, i])
|
||
|
|
||
|
|
||
|
m1.Const_ESS_4 = py.Constraint(m1.h, m1.t, rule=Power_ESS_Disch_limit1,
|
||
|
doc='Discharging power of ESS (Upper limit)')
|
||
|
|
||
|
|
||
|
def Power_ESS_Disch_limit2(m1, h, i):
|
||
|
# return (m1.P_ESS_Disch[h, i]) >= 0.2 * (m1.DischRate_ESS[h] * m1.S_ESS_Disch[h, i])
|
||
|
return (m1.P_ESS_Disch_Home[h, i] + m1.P_ESS_Disch_Grid[h, i]) >= 0
|
||
|
|
||
|
|
||
|
m1.Const_ESS_5 = py.Constraint(m1.h, m1.t, rule=Power_ESS_Disch_limit2,
|
||
|
doc='Discharging power of ESS (Lower limit)')
|
||
|
|
||
|
|
||
|
def Status_ESS(m1, h, i):
|
||
|
return m1.S_ESS_Disch[h, i] + m1.S_ESS_Charge[h, i] <= 1
|
||
|
|
||
|
|
||
|
m1.Const_ESS_6 = py.Constraint(m1.h, m1.t, rule=Status_ESS,
|
||
|
doc='Charging and discharging will not occur at same time')
|
||
|
|
||
|
|
||
|
def Power_ESS_Disch(m1, h, i):
|
||
|
return (m1.P_ESS_Disch[h, i]) * m1.DischEffc_ESS[h] == (m1.P_ESS_Disch_Home[h, i] + m1.P_ESS_Disch_Grid[h, i])
|
||
|
|
||
|
|
||
|
m1.Const_ESS_7 = py.Constraint(m1.h, m1.t, rule=Power_ESS_Disch,
|
||
|
doc='Discharging power of ESS to Home and Grid')
|
||
|
|
||
|
|
||
|
def ESS_SoC_limit1(m1, h, i):
|
||
|
return m1.Energy_ESS[h, i] >= 0.1 * m1.Cap_ESS[h]
|
||
|
|
||
|
|
||
|
m1.Const_ESS_8 = py.Constraint(m1.h, m1.t, rule=ESS_SoC_limit1, doc='Minimum SoC of ESS')
|
||
|
|
||
|
|
||
|
def ESS_SoC_limit2(m1, h, i):
|
||
|
return m1.Energy_ESS[h, i] <= m1.Cap_ESS[h]
|
||
|
|
||
|
|
||
|
m1.Const_ESS_9 = py.Constraint(m1.h, m1.t, rule=ESS_SoC_limit2, doc='Maximum SoC of ESS')
|
||
|
|
||
|
|
||
|
def ESS_final_SoC(m1, h, i):
|
||
|
return m1.Energy_ESS[h, i] >= m1.End_En_ESS[h]
|
||
|
|
||
|
|
||
|
m1.Const_ESS_10 = py.Constraint(m1.h, [m1.t_final.value], rule=ESS_final_SoC,
|
||
|
doc='Final SoC of ESS at departure time')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Electric Water Heater (EWH) Constraints">
|
||
|
def EWH_limit1(m1, h, i):
|
||
|
return m1.tetta_EWH_wat[h, i] <= m1.tetta_up[h]
|
||
|
|
||
|
|
||
|
m1.Const_EWH_1 = py.Constraint(m1.h, m1.t, rule=EWH_limit1, doc='Maximum Limit')
|
||
|
|
||
|
|
||
|
def EWH_limit2(m1, h, i):
|
||
|
return m1.tetta_EWH_wat[h, i] >= m1.tetta_low[h]
|
||
|
|
||
|
|
||
|
m1.Const_EWH_2 = py.Constraint(m1.h, m1.t, rule=EWH_limit2, doc='Minimum Limit')
|
||
|
|
||
|
|
||
|
def EWH_power(m1, h, i):
|
||
|
return m1.P_EWH[h, i] == m1.Q[h] * m1.S_EWH[h, i]
|
||
|
|
||
|
|
||
|
m1.Const_EWH_3 = py.Constraint(m1.h, m1.t, rule=EWH_power, doc='Electric water heater power')
|
||
|
|
||
|
|
||
|
def EWH_temp_1(m1, h, i):
|
||
|
return m1.tetta_EWH_wat[h, i] == m1.tetta_amb[h, i] + m1.Q[h] * m1.S_EWH[h, i] * m1.R[
|
||
|
h] * m1.time_d - (
|
||
|
(m1.M[h] - m1.water_use[h, i]) / m1.M[h]) * (
|
||
|
m1.tetta_amb_int[h] - m1.tetta_EWH_int[h]) * py.exp(-m1.time_d / (m1.R[h] * m1.C[h]))
|
||
|
|
||
|
|
||
|
m1.Const_EWH_4 = py.Constraint(m1.h, [m1.t_first.value], rule=EWH_temp_1, doc='EWH m1')
|
||
|
|
||
|
|
||
|
def EWH_temp_2(m1, h, i):
|
||
|
return m1.tetta_EWH_wat[h, i] == m1.tetta_amb[h, i] + m1.Q[h] * m1.S_EWH[h, i] * m1.R[
|
||
|
h] * m1.time_d - (
|
||
|
(m1.M[h] - m1.water_use[h, i]) / m1.M[h]) * (
|
||
|
m1.tetta_amb[h, i] - m1.tetta_EWH_wat[h, i - 1]) * py.exp(
|
||
|
-m1.time_d / (m1.R[h] * m1.C[h]))
|
||
|
|
||
|
|
||
|
m1.Const_EWH_5 = py.Constraint(m1.h, m1.t_second, rule=EWH_temp_2, doc='EWH m1')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="PV Constraints"> # need to include the PV equation of omer paper
|
||
|
|
||
|
def PV_production(m1, h, i):
|
||
|
return m1.PV_Grid[h, i] + m1.PV_Home[h, i] + m1.PV_Battery[h, i] == m1.PV[h, i]
|
||
|
# return m1.PV_Grid[h, i] + m1.PV_Battery[h, i] == m1.PV[h, i]
|
||
|
|
||
|
|
||
|
m1.Const_PV = py.Constraint(m1.h, m1.t, rule=PV_production, doc='PV Production')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
|
||
|
# <editor-fold desc="objective functions">
|
||
|
|
||
|
|
||
|
def objective_rule(m1):
|
||
|
return sum(
|
||
|
m1.P_Buy_Grid[y, i] * m1.Buy_price[i] - m1.P_Sell_Grid[y, i] * m1.Sell_price[i] for i in
|
||
|
m1.t for y in m1.h) * m1.time_d.value
|
||
|
|
||
|
# return sum(
|
||
|
# (m1.P_Buy_Grid[y, i] * m1.Buy_price[i]) * m1.S_P_buy[y, i] - (m1.P_Sell_Grid[y, i] * m1.Sell_price[i])* m1.S_P_sell[y, i] for i in
|
||
|
# m1.t for y in m1.h) * m1.time_d.value
|
||
|
|
||
|
|
||
|
m1.obj1 = py.Objective(rule=objective_rule, sense=py.minimize, doc='Definition of objective function')
|
||
|
# <editor-fold desc="Solver">
|
||
|
# m1.write('m12.lp', io_options={'symbolic_solver_labels': True})
|
||
|
opt = py.SolverFactory('gurobi')
|
||
|
# opt = py.SolverFactory('mindtpy').solve(m1, mip_solver='cplex', nlp_solver='ipopt', tee=True)
|
||
|
# opt = py.SolverFactory('couenne ')
|
||
|
opt.options["mipgap"] = 0.03 # 0.155 for load and 0.8 for other load
|
||
|
result = opt.solve(m1, tee=True)
|
||
|
# result = opt.solve(m1)
|
||
|
|
||
|
# this is for calling the neos
|
||
|
# os.environ['NEOS_EMAIL'] = 'sadamengr15@gmail.com'
|
||
|
# opt = py.SolverManagerFactory('neos')
|
||
|
# result = opt.solve(m1, opt='cplex')
|
||
|
|
||
|
# m1.pprint()
|
||
|
print(result)
|
||
|
for y in m1.h:
|
||
|
m1.cost[y] = py.value(sum(m1.P_Buy_Grid[y, i] * m1.Buy_price[i] - m1.P_Sell_Grid[y, i] *
|
||
|
m1.Sell_price[i] for i in m1.t) * m1.time_d.value)
|
||
|
for h in m1.h:
|
||
|
print(
|
||
|
f"Before Cost Objective of Home {h} : " f"{py.value((sum((m1.P_Buy_Grid[h, i] * m1.Buy_price[i]) - (m1.P_Sell_Grid[h, i] * m1.Sell_price[i]) for i in m1.t) * m1.time_d.value))}")
|
||
|
# </editor-fold>
|
||
|
|
||
|
# total power on the transformer:
|
||
|
for i in m1.t:
|
||
|
m1.Power_Total_power[i] = round(m1.Power_P_T[i].value + m1.Power_T_P[i].value, 3)
|
||
|
|
||
|
# print('Sum of all homes objectives = ', py.value(m1.obj1))
|
||
|
# for h in m1.h:
|
||
|
# print(
|
||
|
# f"Before Cost Objective of Hom {h} : " f"{py.value((sum((m1.P_Buy_Grid[h, i] * m1.Buy_price[i]) - (m1.P_Sell_Grid[h, i] * m1.Sell_price[i]) for i in m1.t) * m1.time_d))}")
|
||
|
# < / editor - fold >
|
||
|
# for h in m1.h:
|
||
|
# print(
|
||
|
# f"Energy of Home {h} : " f"{py.value((sum((m1.P_Buy_Grid[h, i] * m1.time_d) for i in m1.t)))}")
|
||
|
# for h in m1.h:
|
||
|
# m1.E_Buy_Grid[h] = ((sum((m1.P_Buy_Grid[h, i] * m1.time_d) for i in m1.t) ))
|
||
|
print('power export ', py.value(sum(m1.P_Sell_Grid[y, i] for i in m1.t for y in m1.h)))
|
||
|
print('power import ', py.value(sum(m1.P_Buy_Grid[y, i] for i in m1.t for y in m1.h)))
|
||
|
|
||
|
# # # <editor-fold desc="Ploting">
|
||
|
##
|
||
|
# # # Plotting for automatic ploting
|
||
|
time = [i for i in m1.t]
|
||
|
|
||
|
pl1_pbuy = [[] for j in m1.h]
|
||
|
pl1_pbuy_extra = [[] for j in m1.h]
|
||
|
pl1_pbuy_extra_1 = [[] for j in m1.h]
|
||
|
pl1_psell_extra = [[] for j in m1.h]
|
||
|
pl1_psell_extra_1 = [[] for j in m1.h]
|
||
|
pl_Buy_Price_p2p = [[] for j in m1.h]
|
||
|
pl_Sell_Price_p2p = [[] for j in m1.h]
|
||
|
pl1_pbuy_peer = [[] for j in m1.h]
|
||
|
pl1_psell = [[] for j in m1.h]
|
||
|
pl1_psell_peer = [[] for j in m1.h]
|
||
|
|
||
|
PEVCharge = [[] for j in m1.h]
|
||
|
PESSChargeGrid = [[] for j in m1.h]
|
||
|
PEVDischHome = [[] for j in m1.h]
|
||
|
PESSDischHome = [[] for j in m1.h]
|
||
|
PESSCharge = [[] for j in m1.h]
|
||
|
SEVCharge = [[] for j in m1.h]
|
||
|
SEVDisch = [[] for j in m1.h]
|
||
|
PEVDischGrid = [[] for j in m1.h]
|
||
|
PESSDischGrid = [[] for j in m1.h]
|
||
|
PVHome = [[] for j in m1.h]
|
||
|
PVGrid = [[] for j in m1.h]
|
||
|
PVBattery = [[] for j in m1.h]
|
||
|
PPV = [[] for j in m1.h]
|
||
|
baseload = [[] for j in m1.h]
|
||
|
EnergyESS = [[] for j in m1.h]
|
||
|
EnergyEV = [[] for j in m1.h]
|
||
|
tettaEWHwat = [[] for j in m1.h]
|
||
|
SEWH = [[] for j in m1.h]
|
||
|
PEWH = [[] for j in m1.h]
|
||
|
Buyprice = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
HST = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
DHST = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
DTOT = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
Transformer_Limit = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
Buyprice_peer = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
Power_P_T = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
Power_T_P = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
Sellprice = []
|
||
|
combine = []
|
||
|
cost_1 = [[] for j in m1.h]
|
||
|
Sellprice_peer = []
|
||
|
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_Buy_Grid.items():
|
||
|
if k[0] == j:
|
||
|
pl1_pbuy[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_Buy_extra.items():
|
||
|
if k[0] == j:
|
||
|
pl1_pbuy_extra[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_Buy_extra_1.items():
|
||
|
if k[0] == j:
|
||
|
pl1_pbuy_extra_1[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_sell_extra.items():
|
||
|
if k[0] == j:
|
||
|
pl1_psell_extra[j - 1].append(py.value(v))
|
||
|
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_sell_extra_1.items():
|
||
|
if k[0] == j:
|
||
|
pl1_psell_extra_1[j - 1].append(py.value(v))
|
||
|
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.Buy_Price_p2p.items():
|
||
|
if k[0] == j:
|
||
|
pl_Buy_Price_p2p[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.Sell_Price_p2p.items():
|
||
|
if k[0] == j:
|
||
|
pl_Sell_Price_p2p[j - 1].append(py.value(v))
|
||
|
|
||
|
# for j in m1.h:
|
||
|
# for k, v in m1.P_Sell_Peer.items():
|
||
|
# if k[0] == j:
|
||
|
# pl1_psell_peer[j - 1].append(py.value(v))
|
||
|
# for j in m1.h:
|
||
|
# for k, v in m1.P_Buy_Peer.items():
|
||
|
# if k[0] == j:
|
||
|
# pl1_pbuy_peer[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_Sell_Grid.items():
|
||
|
if k[0] == j:
|
||
|
pl1_psell[j - 1].append(py.value(v))
|
||
|
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_EV_Charge.items():
|
||
|
if k[0] == j:
|
||
|
PEVCharge[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_EV_Disch_Home.items():
|
||
|
if k[0] == j:
|
||
|
PEVDischHome[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_EV_Disch_Grid.items():
|
||
|
if k[0] == j:
|
||
|
PEVDischGrid[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_ESS_Charge_Grid.items():
|
||
|
if k[0] == j:
|
||
|
PESSChargeGrid[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_ESS_Disch_Home.items():
|
||
|
if k[0] == j:
|
||
|
PESSDischHome[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_ESS_Charge.items():
|
||
|
if k[0] == j:
|
||
|
PESSCharge[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_ESS_Disch_Grid.items():
|
||
|
if k[0] == j:
|
||
|
PESSDischGrid[j - 1].append(py.value(v))
|
||
|
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.PV_Home.items():
|
||
|
if k[0] == j:
|
||
|
PVHome[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.PV_Grid.items():
|
||
|
if k[0] == j:
|
||
|
PVGrid[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.PV_Battery.items():
|
||
|
if k[0] == j:
|
||
|
PVBattery[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.PV.items():
|
||
|
if k[0] == j:
|
||
|
PPV[j - 1].append(py.value(v))
|
||
|
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.base_load.items():
|
||
|
if k[0] == j:
|
||
|
baseload[j - 1].append(py.value(v))
|
||
|
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.Energy_ESS.items():
|
||
|
if k[0] == j:
|
||
|
EnergyESS[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.Energy_EV.items():
|
||
|
if k[0] == j:
|
||
|
EnergyEV[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.P_EWH.items():
|
||
|
if k[0] == j:
|
||
|
PEWH[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.tetta_EWH_wat.items():
|
||
|
if k[0] == j:
|
||
|
tettaEWHwat[j - 1].append(py.value(v))
|
||
|
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.S_EWH.items():
|
||
|
if k[0] == j:
|
||
|
SEWH[j - 1].append(py.value(v))
|
||
|
for j in m1.h:
|
||
|
for k, v in m1.cost.items():
|
||
|
if k == j:
|
||
|
cost_1[j - 1].append(py.value(v))
|
||
|
# for i in m1.Buy_price_Peer:
|
||
|
# Buyprice_peer.append(py.value(m1.Buy_price_Peer[i]))
|
||
|
for i in m1.Buy_price:
|
||
|
Buyprice.append(py.value(m1.Buy_price[i]))
|
||
|
|
||
|
for i in m1.Sell_price:
|
||
|
Sellprice.append(py.value(m1.Sell_price[i]))
|
||
|
# for i in m1.Sell_price_Peer:
|
||
|
# Sellprice_peer.append(py.value(m1.Sell_price_Peer[i]))
|
||
|
|
||
|
# # alone
|
||
|
# for k in range(len(m1.h)):
|
||
|
# fig, ax = plt.subplots(5, 3, figsize=(10, 10))
|
||
|
# ax[0, 0].bar(time, pbuy[k], label='Buying power')
|
||
|
# ax[0, 0].bar(time, psell[k], label='Selling power')
|
||
|
# ax[0, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[0, 1].bar(time, baseload[k], label='Base load', color='r')
|
||
|
# ax[0, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[0, 2].plot(time, Buyprice, label='Buyprice')
|
||
|
# ax[0, 2].plot(time, Sellprice, label='Sellprice')
|
||
|
# ax[0, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[1, 0].bar(time, PEWH[k], label='EWH Power', color='r')
|
||
|
# ax[1, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[1, 1].plot(time, tettaEWHwat[k], label='EWH Temp')
|
||
|
# ax[1, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[1, 2].bar(time, SEWH[k], label='Status of EWH')
|
||
|
# ax[1, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[2, 0].bar(time, PESSChargeGrid[k], label='ESS Charging from Grid', color='r')
|
||
|
# ax[2, 0].bar(time, PVBattery[k], label='ESS Charging from PV', color='g')
|
||
|
# ax[2, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[2, 1].bar(time, PESSDischHome[k], label='ESS Disch to home')
|
||
|
# ax[2, 1].bar(time, PESSDischGrid[k], label='ESS Disch to grid')
|
||
|
# ax[2, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[2, 2].plot(time, EnergyESS[k], label='Energy of ESS', color='g')
|
||
|
# ax[2, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[3, 0].bar(time, PVHome[k], label='PV to Home')
|
||
|
# ax[3, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[3, 1].bar(time, PVGrid[k], label='PV to Grid')
|
||
|
# ax[3, 1].bar(time, PVBattery[k], label='PV to battery')
|
||
|
# ax[3, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[3, 2].bar(time, PVHome[k], label='PV to Home')
|
||
|
# ax[3, 2].bar(time, PVGrid[k], label='PV to Grid')
|
||
|
# ax[3, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[4, 0].bar(time, PEVCharge[k], label='EV Charging power', color='r')
|
||
|
# ax[4, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[4, 1].bar(time, PEVDischHome[k], label='EV Disch to home')
|
||
|
# ax[4, 1].bar(time, PEVDischGrid[k], label='EV Disch to grid')
|
||
|
# ax[4, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[4, 2].plot(time, EnergyEV[k], label='Energy of EV', color='r')
|
||
|
# ax[4, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# ax[4, 0].set_xlabel('Time (step)')
|
||
|
# ax[4, 1].set_xlabel('Time (step)')
|
||
|
# ax[4, 2].set_xlabel('Time (step)')
|
||
|
# ax[0, 2].set_ylabel('price ($/W/h)')
|
||
|
# ax[0, 0].set_ylabel('Power (kW)')
|
||
|
# ax[1, 0].set_ylabel('Power (kW)')
|
||
|
# ax[2, 0].set_ylabel('Power (kW)')
|
||
|
# ax[3, 0].set_ylabel('Power (kW)')
|
||
|
# ax[4, 0].set_ylabel('Power (kW)')
|
||
|
# ax[2, 2].set_ylabel('Energy (kWh)')
|
||
|
# ax[4, 2].set_ylabel('Energy (kWh)')
|
||
|
# plt.suptitle(
|
||
|
# f" Home {k + 1} : " f" Cost is : {py.value((sum((m1.P_Buy_Grid[k + 1, i] * m1.Buy_price[i]) - (m1.P_Sell_Grid[k + 1, i] * m1.Sell_price[i]) for i in m1.t) * m1.time_d))}")
|
||
|
# pl.tight_layout()
|
||
|
# plt.show()
|
||
|
for k in m1.t:
|
||
|
m1.combine[k] = (PEVCharge[1][k - 1] + PEWH[1][k - 1] + baseload[1][k - 1] + PESSChargeGrid[1][k - 1])
|
||
|
# # code for if you want to plot all the home togathoer
|
||
|
for i in m1.combine:
|
||
|
combine.append(py.value(m1.combine[i]))
|
||
|
# fig, ax = plt.subplots(6, len(m1.h), constrained_layout=True, sharex='col',sharey='row') # , sharey='row'
|
||
|
# gs = fig.add_gridspec(hspace=0, wspace=0)
|
||
|
# for k in range(len(m1.h)):
|
||
|
# ax[0, k].bar(time, pl1_pbuy[k], color="blue")
|
||
|
# ax[1, k].bar(time, pl1_psell[k], color="black")
|
||
|
# # ax[2, k].bar(time, pl1_pbuy_peer[k], color="green")
|
||
|
# # ax[3, k].bar(time, pl1_psell_peer[k], color="red")
|
||
|
# # ax[4, 0].plot(time, Transformer_Limit, color="red")
|
||
|
# # ax[4, 0].bar(time, Power_P_T, color="blue")
|
||
|
# # ax[4, 0].legend(['Tranformer_limit'])
|
||
|
# # ax[4, 1].plot(time, Transformer_Limit, color="red")
|
||
|
# # ax[4, 1].bar(time, Power_T_P, color="blue")
|
||
|
# # ax[4, 1].legend(['Power_T_P'])
|
||
|
# # # ax[0, k].plot(time, baseload[k], label='Base load')
|
||
|
#
|
||
|
# ax[2, k].plot(time, Buyprice, label='Buyprice_G', color="blue")
|
||
|
# # ax[4, 2].plot(time, Buyprice_peer, label='Buyprice_P', color="red")
|
||
|
# ax[2, k].plot(time, Sellprice, label='Sellprice_G', color="green")
|
||
|
# # ax[4, 2].plot(time, Sellprice_peer, label='Sellprice_P', color="black")
|
||
|
# ax[2, 0].legend(['Buyprice_G', 'Sellprice_G'],loc='upper left', fontsize='xx-small')
|
||
|
# # ax[4, 3].plot(time, HST, label='HST', color="red")
|
||
|
# # ax[4, 3].plot(time, DHST, label='DHST', color="blue")
|
||
|
# # ax[4, 3].plot(time, DTOT, label='DTOT', color="green")
|
||
|
# # ax[4, 3].legend(['HST', 'DHST', 'DTOT'])
|
||
|
#
|
||
|
# # ax[9, k].plot(time, PEWH[k], label='EWH Power')
|
||
|
# ax[3, k].plot(time, tettaEWHwat[k], label='EWH Temp')
|
||
|
# ax[3, 0].legend(ncol=3,loc='upper left', fontsize='xx-small')
|
||
|
# # ax[5, 0].legend(ncol=3,loc='best', fontsize='xx-small')
|
||
|
# # ax[3, k].plot(time, PESSCharge[k], label='ESS Charging power')
|
||
|
# # ax[3, k].plot(time, PESSDischHome[k], label='ESS Disch to home')
|
||
|
# # ax[3, k].plot(time, PESSDischGrid[k], label='ESS Disch to grid')
|
||
|
# # ax[3, 0].legend(ncol=3,loc='best', fontsize='xx-small')
|
||
|
# # ax[4, k].plot(time, PVHome[k], label='PV to Home')
|
||
|
# # ax[4, k].plot(time, PVGrid[k], label='PV to Grid')
|
||
|
# # ax[4, k].plot(time, PPV[k], label='All PV Production')
|
||
|
# # ax[4, 0].legend(ncol=3,loc='best', fontsize='xx-small')
|
||
|
# # ax[5, k].plot(time, PEVCharge[k], label='EV Charging power')
|
||
|
# # ax[5, k].plot(time, PEVDischHome[k], label='EV Disch to home')
|
||
|
# # ax[5, k].plot(time, PEVDischGrid[k], label='EV Disch to grid')
|
||
|
# # ax[5, 0].legend(ncol=3,loc='best', fontsize='xx-small')
|
||
|
# ax[4, k].plot(time, EnergyESS[k], label='Energy of ESS')
|
||
|
# ax[4, 0].legend(ncol=3,loc='upper left', fontsize='xx-small')
|
||
|
# ax[5, k].plot(time, EnergyEV[k], label='Energy of EV')
|
||
|
# ax[5, 0].legend(ncol=3,loc='upper left', fontsize='xx-small')
|
||
|
# ax[5, k].set_xlabel('Time (step)')
|
||
|
# # # ax[1, k].set_ylabel('Electricity price ($/W/h)')
|
||
|
# ax[0, 0].set_ylabel('Dem', fontsize='x-small')
|
||
|
# ax[1, 0].set_ylabel('Sel', fontsize='x-small')
|
||
|
# ax[0, k].set_title(f"Home {k + 1}", fontsize='xx-small')
|
||
|
# fig.suptitle('Home Profile before limit', fontsize=16)
|
||
|
# # #
|
||
|
# # pl.tight_layout()
|
||
|
#
|
||
|
# plt.show()
|
||
|
|
||
|
# # if you want to plot each home alone in the figure
|
||
|
# # for k in range(len(m1.h)):
|
||
|
# # fig, ax = plt.subplots(2)
|
||
|
# # ax[0].plot(time, pbuy[k], 'b-', label='Base')
|
||
|
# # ax[0].plot(time, psell[k], 'g--', label='Controllable Load')
|
||
|
# # ax[1].plot(time, Buyprice, label='buying price')
|
||
|
# # ax[1].plot(time, Sellprice, 'y-+', label='selling price')
|
||
|
# # ax[1].legend(loc='upper left')
|
||
|
# # ax[1].set_xlabel('Time (hour)')
|
||
|
# # ax[1].set_ylabel('Electricity price ($/W/h)')
|
||
|
# # plt.suptitle(f"Home {k+1}")
|
||
|
# # pl.tight_layout()
|
||
|
# # plt.show()
|
||
|
#
|
||
|
# print('this is the end of code')
|
||
|
# # # # </editor-fold>
|
||
|
|
||
|
# </editor-fold>
|
||
|
# #</editor-fold>
|
||
|
|
||
|
#
|
||
|
# # <editor-fold desc="Ploting">
|
||
|
|
||
|
#
|
||
|
# # # # <editor-fold desc="Ploting">
|
||
|
# # #
|
||
|
# # # # Plotting for automatic ploting
|
||
|
# # time = [i for i in m1.t]
|
||
|
# pl2_pbuy = [[] for j in m1.h]
|
||
|
# # pl2_psell = [[] for j in m1.h]
|
||
|
# # # # PEVCharge = [[] for j in m1.h]
|
||
|
# # # # PESSChargeGrid = [[] for j in m1.h]
|
||
|
# # # # PEVDischHome = [[] for j in m1.h]
|
||
|
# # # # PESSDischHome = [[] for j in m1.h]
|
||
|
# # # # SEVCharge = [[] for j in m1.h]
|
||
|
# # # # SEVDisch = [[] for j in m1.h]
|
||
|
# # # # PEVDischGrid = [[] for j in m1.h]
|
||
|
# # # # PESSDischGrid = [[] for j in m1.h]
|
||
|
# # # # PVHome = [[] for j in m1.h]
|
||
|
# # # # PVGrid = [[] for j in m1.h]
|
||
|
# # # # PVBattery = [[] for j in m1.h]
|
||
|
# # # # PPV = [[] for j in m1.h]
|
||
|
# # # # baseload = [[] for j in m1.h]
|
||
|
# # # # EnergyESS = [[] for j in m1.h]
|
||
|
# # # # EnergyEV = [[] for j in m1.h]
|
||
|
# # # # tettaEWHwat = [[] for j in m1.h]
|
||
|
# # # # SEWH = [[] for j in m1.h]
|
||
|
# # # # PEWH = [[] for j in m1.h]
|
||
|
# # # # Buyprice = [] # cost = [i for i in m1.c.values()] <<< check this
|
||
|
# # # # Sellprice = []
|
||
|
#
|
||
|
# for j in m1.h:
|
||
|
# for k, v in m1.P_Buy_Grid.items():
|
||
|
# if k[0] == j:
|
||
|
# pl2_pbuy[j - 1].append(py.value(v))
|
||
|
# # for j in m1.h:
|
||
|
# # for k, v in m1.P_Sell_Grid.items():
|
||
|
# # if k[0] == j:
|
||
|
# # pl2_psell[j - 1].append(py.value(v))
|
||
|
# # # #
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.P_EV_Charge.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PEVCharge[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.P_EV_Disch_Home.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PEVDischHome[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.P_EV_Disch_Grid.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PEVDischGrid[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.P_ESS_Charge_Grid.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PESSChargeGrid[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.P_ESS_Disch_Home.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PESSDischHome[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.P_ESS_Disch_Grid.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PESSDischGrid[j - 1].append(py.value(v))
|
||
|
# # # #
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.PV_Home.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PVHome[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.PV_Grid.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PVGrid[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.PV_Battery.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PVBattery[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.PV.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PPV[j - 1].append(py.value(v))
|
||
|
# # # #
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.base_load.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # baseload[j - 1].append(py.value(v))
|
||
|
# # # #
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.Energy_ESS.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # EnergyESS[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.Energy_EV.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # EnergyEV[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.P_EWH.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # PEWH[j - 1].append(py.value(v))
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.tetta_EWH_wat.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # tettaEWHwat[j - 1].append(py.value(v))
|
||
|
# # # #
|
||
|
# # # # for j in m1.h:
|
||
|
# # # # for k, v in m1.S_EWH.items():
|
||
|
# # # # if k[0] == j:
|
||
|
# # # # SEWH[j - 1].append(py.value(v))
|
||
|
# # # #
|
||
|
# # # # for i in m1.Buy_price:
|
||
|
# # # # Buyprice.append(value(m1.Buy_price[i]))
|
||
|
# # # #
|
||
|
# # # # for i in m1.Sell_price:
|
||
|
# # # # Sellprice.append(value(m1.Sell_price[i]))
|
||
|
# # # #
|
||
|
# # # # # alone
|
||
|
# # # for k in range(len(m1.h)):
|
||
|
# # # fig, ax = plt.subplots(5, 2, figsize=(10, 10))
|
||
|
# # # ax[0, 0].bar(time, pl2_pbuy[k], label='Buying power-obj2')
|
||
|
# # # ax[0, 0].bar(time, pl1_pbuy[k], label='Buying power-obj1')
|
||
|
# # # ax[1, 0].bar(time, pl2_psell[k], label='Selling power-obj2')
|
||
|
# # # ax[1, 0].bar(time, pl1_psell[k], label='Selling power-obj1')
|
||
|
# # # ax[0, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # ax[1, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[0, 1].bar(time, baseload[k], label='Base load', color='r')
|
||
|
# # # # # ax[0, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[0, 2].plot(time, Buyprice, label='Buyprice')
|
||
|
# # # # # ax[0, 2].plot(time, Sellprice, label='Sellprice')
|
||
|
# # # # # ax[0, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[1, 0].bar(time, PEWH[k], label='EWH Power', color='r')
|
||
|
# # # # # ax[1, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[1, 1].plot(time, tettaEWHwat[k], label='EWH Temp')
|
||
|
# # # # # ax[1, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[1, 2].bar(time, SEWH[k], label='Status of EWH')
|
||
|
# # # # # ax[1, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[2, 0].bar(time, PESSChargeGrid[k], label='ESS Charging from Grid', color='r')
|
||
|
# # # # # ax[2, 0].bar(time, PVBattery[k], label='ESS Charging from PV', color='g')
|
||
|
# # # # # ax[2, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[2, 1].bar(time, PESSDischHome[k], label='ESS Disch to home')
|
||
|
# # # # # ax[2, 1].bar(time, PESSDischGrid[k], label='ESS Disch to grid')
|
||
|
# # # # # ax[2, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[2, 2].plot(time, EnergyESS[k], label='Energy of ESS', color='g')
|
||
|
# # # # # ax[2, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[3, 0].bar(time, PVHome[k], label='PV to Home')
|
||
|
# # # # # ax[3, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[3, 1].bar(time, PVGrid[k], label='PV to Grid')
|
||
|
# # # # # ax[3, 1].bar(time, PVBattery[k], label='PV to battery')
|
||
|
# # # # # ax[3, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[3, 2].bar(time, PVHome[k], label='PV to Home')
|
||
|
# # # # # ax[3, 2].bar(time, PVGrid[k], label='PV to Grid')
|
||
|
# # # # # ax[3, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[4, 0].bar(time, PEVCharge[k], label='EV Charging power', color='r')
|
||
|
# # # # # ax[4, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[4, 1].bar(time, PEVDischHome[k], label='EV Disch to home')
|
||
|
# # # # # ax[4, 1].bar(time, PEVDischGrid[k], label='EV Disch to grid')
|
||
|
# # # # # ax[4, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[4, 2].plot(time, EnergyEV[k], label='Energy of EV', color='r')
|
||
|
# # # # # ax[4, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # # ax[4, 0].set_xlabel('Time (step)')
|
||
|
# # # # # ax[4, 1].set_xlabel('Time (step)')
|
||
|
# # # # # ax[4, 2].set_xlabel('Time (step)')
|
||
|
# # # # # ax[0, 2].set_ylabel('price ($/W/h)')
|
||
|
# # # # # ax[0, 0].set_ylabel('Power1 (kW)')
|
||
|
# # # # # ax[1, 0].set_ylabel('Power (kW)')
|
||
|
# # # # # ax[2, 0].set_ylabel('Power (kW)')
|
||
|
# # # # # ax[3, 0].set_ylabel('Power (kW)')
|
||
|
# # # # # ax[4, 0].set_ylabel('Power (kW)')
|
||
|
# # # # # ax[2, 2].set_ylabel('Energy (kWh)')
|
||
|
# # # # # ax[4, 2].set_ylabel('Energy (kWh)')
|
||
|
# # # plt.suptitle(
|
||
|
# # # f" Home {k + 1} : " f" Cost is : {py.value((sum(m1.P_Buy_Grid[k + 1, i] for i in m1.t) * m1.time_d))}")
|
||
|
# # # pl.tight_layout()
|
||
|
# # # plt.show()
|
||
|
# # # #
|
||
|
# # # # # code for if you want to plot all the home togathoer
|
||
|
# # fig, ax = plt.subplots(4, len(m1.h), figsize=(10, 10), constrained_layout=True, sharex='col', sharey='row')
|
||
|
# # gs = fig.add_gridspec(hspace=0, wspace=0)
|
||
|
# # for k in range(len(m1.h)):
|
||
|
# # ax[0, k].bar(time, pl1_pbuy[k], color='green')
|
||
|
# # ax[1, k].bar(time, pl2_pbuy[k], color='blue')
|
||
|
# # ax[2, k].bar(time, pl1_psell[k], color='green')
|
||
|
# # ax[3, k].bar(time, pl2_psell[k], label='blue')
|
||
|
# # # ax[0, k].legend(loc='best')
|
||
|
# # # ax[1, k].legend(loc='best')
|
||
|
# # # ax[2, k].legend(loc='best')
|
||
|
# # # ax[3, k].legend(loc='best')
|
||
|
# # # # ax[0, k].plot(time, baseload[k], label='Base load')
|
||
|
# # # # # ax[1, k].plot(time, Buyprice, label='Buyprice')
|
||
|
# # # # ax[1, k].legend(loc='best')
|
||
|
# # # # # ax[1, k].plot(time, Sellprice, label='Sellprice')
|
||
|
# # # # # ax[1, k].legend(loc='best')
|
||
|
# # # # # ax[2, k].plot(time, PEWH[k], label='EWH Power')
|
||
|
# # # # # ax[2, k].plot(time, tettaEWHwat[k], label='EWH Temp')
|
||
|
# # # # # ax[2, k].legend(loc='best')
|
||
|
# # # # # ax[3, k].plot(time, PESSCharge[k],label='ESS Charging power')
|
||
|
# # # # # ax[3, k].plot(time, PESSDischHome[k], label='ESS Disch to home')
|
||
|
# # # # # ax[3, k].plot(time, PESSDischGrid[k], label='ESS Disch to grid')
|
||
|
# # # # # ax[3, k].legend(loc='best')
|
||
|
# # # # # ax[4, k].plot(time, PVHome[k], label='PV to Home')
|
||
|
# # # # # ax[4, k].plot(time, PVGrid[k], label='PV to Grid')
|
||
|
# # # # # ax[4, k].plot(time, PPV[k], label='All PV Production')
|
||
|
# # # # # ax[4, k].legend(loc='best')
|
||
|
# # # # # ax[5, k].plot(time, PEVCharge[k], label='EV Charging power')
|
||
|
# # # # # ax[5, k].plot(time, PEVDischHome[k], label='EV Disch to home')
|
||
|
# # # # # ax[5, k].plot(time, PEVDischGrid[k], label='EV Disch to grid')
|
||
|
# # # # # ax[5, k].legend(loc='best')
|
||
|
# # # # # ax[6, k].plot(time, EnergyESS[k], label='Energy of ESS')
|
||
|
# # # # # ax[6, k].legend(loc='best')
|
||
|
# # # # # ax[7, k].plot(time, EnergyEV[k], label='Energy of EV')
|
||
|
# # # # # ax[7, k].legend(loc='best')
|
||
|
# # # # # ax[4, k].set_xlabel('Time (step)')
|
||
|
# # # # # ax[1, k].set_ylabel('Electricity price ($/W/h)')
|
||
|
# # ax[0, 0].set_ylabel(
|
||
|
# # 'Buy-Power1 (kW)') # showing only on the side of the graph if you want to put it on all need to use [0,k]
|
||
|
# # ax[1, 0].set_ylabel('Buy-Power2 (kW)')
|
||
|
# # ax[2, 0].set_ylabel('Sell-Power1 (kW)')
|
||
|
# # ax[3, 0].set_ylabel('Sell-Power2 (kW)')
|
||
|
# #
|
||
|
# # ax[0, k].set_title(f"Home {k + 1}", fontsize='xx-small')
|
||
|
# # fig.suptitle('Cost and Energy minimization', fontsize=16)
|
||
|
# #
|
||
|
# # plt.show()
|
||
|
# # #
|
||
|
# # # # if you want to plot each home alone in the figure
|
||
|
# # # # for k in range(len(m1.h)):
|
||
|
# # # # fig, ax = plt.subplots(2)
|
||
|
# # # # ax[0].plot(time, pbuy[k], 'b-', label='Base')
|
||
|
# # # # ax[0].plot(time, psell[k], 'g--', label='Controllable Load')
|
||
|
# # # # ax[1].plot(time, Buyprice, label='buying price')
|
||
|
# # # # ax[1].plot(time, Sellprice, 'y-+', label='selling price')
|
||
|
# # # # ax[1].legend(loc='upper left')
|
||
|
# # # # ax[1].set_xlabel('Time (hour)')
|
||
|
# # # # ax[1].set_ylabel('Electricity price ($/W/h)')
|
||
|
# # # # plt.suptitle(f"Home {k+1}")
|
||
|
# # # # pl.tight_layout()
|
||
|
# # # # plt.show()
|
||
|
# # #
|
||
|
# # # print('this is the end of code')
|
||
|
# # # </editor-fold>
|
||
|
#
|
||
|
# writer = pd.ExcelWriter('Result_HEMS.xlsx', engine='xlsxwriter')
|
||
|
# time_excel = pd.DataFrame({'Time_Step': time})
|
||
|
# time_excel.to_excel(writer, sheet_name='Power', startcol=0, index=False, header=True)
|
||
|
# for h in range(len(m1.h)):
|
||
|
# pload_excel = pd.DataFrame({f"Power_Demand {h + 1}": pl1_pbuy[h]})
|
||
|
# pload_excel.to_excel(writer, sheet_name='Power', startcol=h + 1, index=False, header=True)
|
||
|
#
|
||
|
# pload_excel = pd.DataFrame({f"Power_Sell {h + 1}": pl1_psell[h]})
|
||
|
# pload_excel.to_excel(writer, sheet_name='Power', startcol=h + 4, index=False, header=True)
|
||
|
# writer.save()
|
||
|
|
||
|
|
||
|
# fig_HEMS, ax = plt.subplots(5, len(m1.h), constrained_layout=True, sharex='col', sharey='row') # , sharey='row'
|
||
|
# gs = fig_HEMS.add_gridspec(hspace=0, wspace=0)
|
||
|
# for k in range(len(m1.h)):
|
||
|
# ax[0, k].bar(time, pl1_pbuy[k], color="blue")
|
||
|
# ax[1, k].bar(time, pl1_psell[k], color="black")
|
||
|
# # ax[2, k].bar(time, pl1_pbuy_peer[k], color="green")
|
||
|
# # ax[3, k].bar(time, pl1_psell_peer[k], color="red")
|
||
|
# # ax[4, 0].plot(time, Transformer_Limit, color="red")
|
||
|
# # ax[4, 0].bar(time, Power_P_T, color="blue")
|
||
|
# # ax[4, 0].legend(['Tranformer_limit'])
|
||
|
# # ax[4, 1].plot(time, Transformer_Limit, color="red")
|
||
|
# # ax[4, 1].bar(time, Power_T_P, color="blue")
|
||
|
# # ax[4, 1].legend(['Power_T_P'])
|
||
|
# # # ax[0, k].plot(time, baseload[k], label='Base load')
|
||
|
#
|
||
|
# # ax[2, k].plot(time, Buyprice, label='Buyprice_G', color="blue")
|
||
|
# # ax[4, 2].plot(time, Buyprice_peer, label='Buyprice_P', color="red")
|
||
|
# # ax[2, k].plot(time, Sellprice, label='Sellprice_G', color="green")
|
||
|
# # ax[4, 2].plot(time, Sellprice_peer, label='Sellprice_P', color="black")
|
||
|
# # ax[2, 0].legend(['Buyprice_G', 'Sellprice_G'],loc='upper left', fontsize='xx-small')
|
||
|
# # ax[4, 3].plot(time, HST, label='HST', color="red")
|
||
|
# # ax[4, 3].plot(time, DHST, label='DHST', color="blue")
|
||
|
# # ax[4, 3].plot(time, DTOT, label='DTOT', color="green")
|
||
|
# # ax[4, 3].legend(['HST', 'DHST', 'DTOT'])
|
||
|
#
|
||
|
# # ax[9, k].plot(time, PEWH[k], label='EWH Power')
|
||
|
# ax[2, k].plot(time, tettaEWHwat[k], label='EWH Temp')
|
||
|
# ax[2, 0].legend(ncol=3, loc='upper left', fontsize='xx-small')
|
||
|
# # ax[5, 0].legend(ncol=3,loc='best', fontsize='xx-small')
|
||
|
# # ax[3, k].plot(time, PESSCharge[k], label='ESS Charging power')
|
||
|
# # ax[3, k].plot(time, PESSDischHome[k], label='ESS Disch to home')
|
||
|
# # ax[3, k].plot(time, PESSDischGrid[k], label='ESS Disch to grid')
|
||
|
# # ax[3, 0].legend(ncol=3,loc='best', fontsize='xx-small')
|
||
|
# # ax[4, k].plot(time, PVHome[k], label='PV to Home')
|
||
|
# # ax[4, k].plot(time, PVGrid[k], label='PV to Grid')
|
||
|
# # ax[4, k].plot(time, PPV[k], label='All PV Production')
|
||
|
# # ax[4, 0].legend(ncol=3,loc='best', fontsize='xx-small')
|
||
|
# # ax[5, k].plot(time, PEVCharge[k], label='EV Charging power')
|
||
|
# # ax[5, k].plot(time, PEVDischHome[k], label='EV Disch to home')
|
||
|
# # ax[5, k].plot(time, PEVDischGrid[k], label='EV Disch to grid')
|
||
|
# # ax[5, 0].legend(ncol=3,loc='best', fontsize='xx-small')
|
||
|
# ax[3, k].plot(time, EnergyESS[k], label='Energy of ESS')
|
||
|
# ax[3, 0].legend(ncol=3, loc='upper left', fontsize='xx-small')
|
||
|
# ax[4, k].bar(time, EnergyEV[k], label='Energy of EV')
|
||
|
# ax[4, 0].legend(ncol=3, loc='upper left', fontsize='xx-small')
|
||
|
# ax[4, k].set_xlabel('Time (step)')
|
||
|
# # # ax[1, k].set_ylabel('Electricity price ($/W/h)')
|
||
|
# ax[0, 0].set_ylabel('Dem', fontsize='x-small')
|
||
|
# ax[1, 0].set_ylabel('Sel', fontsize='x-small')
|
||
|
# ax[0, k].set_title(f"Home {k + 1}", fontsize='xx-small')
|
||
|
# fig_HEMS.suptitle('Home Profile', fontsize=16)
|
||
|
# # #
|
||
|
# # pl.tight_layout()
|
||
|
#
|
||
|
# plt.show()
|
||
|
|
||
|
# stack plot of power demand of one house :
|
||
|
# Generate some example data
|
||
|
# Plot stacked bar chart
|
||
|
fig_all, axs = plt.subplots(5, 2, figsize=(8, 8),constrained_layout=True, sharex='col', sharey='row')
|
||
|
labels = ('EV_Charging', 'Power_EWH', 'Baseload', 'Battery_Charging')
|
||
|
data_lists = [(PEVCharge[i], PEWH[i], baseload[i], PESSChargeGrid[i]) for i in range(10)]
|
||
|
rc('font', family='Times New Roman')
|
||
|
for i, ax in enumerate(axs.flat):
|
||
|
bottom = np.zeros(len(time))
|
||
|
ax.set_xlim([0, 98])
|
||
|
ax.set_ylim([0, 18])
|
||
|
ax.set_xticks(np.arange(0, 98, 6))
|
||
|
# ax.set_yticks(np.arange(0, 24, 4))
|
||
|
for label, data in zip(labels, data_lists[i]):
|
||
|
ax.bar(time, data, bottom=bottom, label=label)
|
||
|
bottom += data
|
||
|
# if i == 3 or i == 6 or i ==9:
|
||
|
# ax.set_title(f'Consumer-{i+1}', fontname='Times New Roman', fontsize=12)
|
||
|
# else:
|
||
|
ax.set_title(f'Prosumer-{i + 1}', fontname='Times New Roman', fontsize=12)
|
||
|
ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
if i == 0 or i == 2 or i == 4 or i == 6 or i == 8:
|
||
|
ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
ax.set_yticks(np.arange(0, 18, 4))
|
||
|
else:
|
||
|
ax.set_ylabel('')
|
||
|
if i == 8 or i == 9:
|
||
|
ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman', fontsize=12)
|
||
|
ax.tick_params(axis='x', which='both', labelsize=12)
|
||
|
else:
|
||
|
ax.tick_params(axis='x', which='both', bottom=True, labelbottom=False) # modified
|
||
|
ax.set_xticks([]) # added
|
||
|
handles, labels = ax.get_legend_handles_labels()
|
||
|
|
||
|
fig_all.legend(handles, labels, loc='lower center', ncol=4, fontsize=12)
|
||
|
plt.tight_layout()
|
||
|
fig_all.subplots_adjust(bottom=0.129)
|
||
|
plt.show()
|
||
|
#
|
||
|
# fig_h1, ax = plt.subplots(constrained_layout=True,figsize=(8, 4))
|
||
|
# labels_h1 = ('EV_Charging', 'Power_EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list_h1 = [PEVCharge[0], PEWH[0], baseload[0], PESSChargeGrid[0]]
|
||
|
# rc('font', family='Times New Roman')
|
||
|
# bottom_h1 = np.zeros(len(time))
|
||
|
# for label_h1, data_h1 in zip(labels_h1, data_list_h1):
|
||
|
# ax.bar(time, data_h1, bottom=bottom_h1, label=label_h1)
|
||
|
# bottom_h1 += data_h1
|
||
|
# # Set plot labels and title
|
||
|
# ax.set_ylabel('Power demand (kW)',fontname='Times New Roman',fontsize=12)
|
||
|
# # ax.set_title('Power demand (Prosumer-1)', fontname='Times New Roman')
|
||
|
# ax.legend(loc='upper left',fontsize=12,ncol=2)
|
||
|
# ax.set_xlabel('Time (15 min interval)',fontname='Times New Roman',fontsize=12)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
# # Show plot
|
||
|
# plt.show()
|
||
|
#
|
||
|
#
|
||
|
# fig_h2, ax = plt.subplots(constrained_layout=True,figsize=(8, 4))
|
||
|
# labels_h2 = ('EV_Charging', 'Power_EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list_h2 = [PEVCharge[1], PEWH[1], baseload[1], PESSChargeGrid[1]]
|
||
|
# rc('font', family='Times New Roman')
|
||
|
# bottom_h2 = np.zeros(len(time))
|
||
|
# for label_h2, data_h2 in zip(labels_h2, data_list_h2):
|
||
|
# ax.bar(time, data_h2, bottom=bottom_h2, label=label_h2)
|
||
|
# bottom_h2 += data_h2
|
||
|
# # ax.plot(time,combine)
|
||
|
# # Set plot labels and title
|
||
|
# ax.set_ylabel('Power demand (kW)', fontname='Times New Roman',fontsize=12)
|
||
|
# # ax.set_title('Power demand (Prosumer-2)', fontname='Times New Roman',fontsize=12)
|
||
|
# ax.legend(loc='upper left',fontsize=12,ncol=2)
|
||
|
# ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman',fontsize=12)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
# # Show plot
|
||
|
# plt.show()
|
||
|
#
|
||
|
#
|
||
|
# fig_h3, ax = plt.subplots(constrained_layout=True,figsize=(8, 4))
|
||
|
# labels_h3 = ( 'Power_EWH', 'Baseload')
|
||
|
# data_list_h3 = [PEWH[2], baseload[2]]
|
||
|
# rc('font', family='Times New Roman')
|
||
|
# colors = ['orange', 'green']
|
||
|
# bottom_h3 = np.zeros(len(time))
|
||
|
# # for i,(label_h3, data_h3) in zip(labels_h3, data_list_h3):
|
||
|
# # ax.bar(time, data_h3, bottom=bottom_h3, label=label_h3, color=colors[label_h3])
|
||
|
# # bottom_h3 += data_h3
|
||
|
# for i, (label_h3, data_h3) in enumerate(zip(labels_h3, data_list_h3)):
|
||
|
# ax.bar(time, data_h3, bottom=bottom_h3, label=label_h3, color=colors[i])
|
||
|
# bottom_h3 += data_h3
|
||
|
#
|
||
|
# # ax.plot(time,combine)
|
||
|
# # Set plot labels and title
|
||
|
# ax.set_ylabel('Power demand (kW)', fontname='Times New Roman',fontsize=12)
|
||
|
# # ax.set_title('Power demand (Consumer-3)', fontname='Times New Roman',fontsize=12)
|
||
|
# ax.legend(loc='upper left',ncol=2,fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman',fontsize=12)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
# # Show plot
|
||
|
# plt.show()
|
||
|
|
||
|
|
||
|
# stack plot for the selling energy to the grid
|
||
|
# Plot stacked bar chart
|
||
|
|
||
|
fig_all_sell, axs = plt.subplots(5, 2, sharex='col',constrained_layout=True, figsize=(8, 8), sharey='row')
|
||
|
labels = ('EV-Discharging_Grid', 'EV-Discharging_Home', 'PV_Grid', 'PV_Home', 'PV_Battery', 'Battery-Discharging_Grid',
|
||
|
'Battery-Discharging_Home')
|
||
|
data_lists = [(PEVDischGrid[k], PEVDischHome[k], PVGrid[k], PVHome[k], PVBattery[k], PESSDischGrid[k], PESSDischHome[k])
|
||
|
for k in range(10)]
|
||
|
rc('font', family='Times New Roman')
|
||
|
for i, ax in enumerate(axs.flat):
|
||
|
bottom = np.zeros(len(time))
|
||
|
ax.set_xlim([0, 98])
|
||
|
ax.set_ylim([0, 14])
|
||
|
ax.set_xticks(np.arange(0, 98, 6))
|
||
|
# ax.set_xticks(
|
||
|
# ['9am', '10am', '11am', '12pm', '1pm', '2pm', '3pm', '4pm', '5pm', '6pm', '7pm', '8pm', '9pm', '10pm', '11pm',
|
||
|
# '12am', '1am','2am', '3am', '4am', '5am', '6am', '7am', '8am'])
|
||
|
|
||
|
for label, data in zip(labels, data_lists[i]):
|
||
|
ax.bar(time, data, bottom=bottom, label=label)
|
||
|
bottom += data
|
||
|
# if i == 3 or i == 6 or i ==9:
|
||
|
# ax.set_title(f'Consumer-{i+1}', fontname='Times New Roman', fontsize=12)
|
||
|
# else:
|
||
|
ax.set_title(f'Prosumer-{i + 1}', fontname='Times New Roman', fontsize=12)
|
||
|
ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
if i == 0 or i == 2 or i == 4 or i == 6 or i == 8:
|
||
|
ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
ax.set_yticks(np.arange(0, 14, 2))
|
||
|
else:
|
||
|
ax.set_ylabel('')
|
||
|
if i == 8 or i == 9:
|
||
|
ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman', fontsize=12)
|
||
|
ax.tick_params(axis='x', which='both', labelsize=12)
|
||
|
else:
|
||
|
ax.tick_params(axis='x', which='both', bottom=True, labelbottom=False) # modified
|
||
|
ax.set_xticks([]) # added
|
||
|
handles, labels = ax.get_legend_handles_labels()
|
||
|
fig_all_sell.legend(handles, labels, loc='lower center', ncol=4, fontsize=10)
|
||
|
plt.tight_layout()
|
||
|
fig_all_sell.subplots_adjust(bottom=0.14)
|
||
|
plt.show()
|
||
|
|
||
|
# fig2_h1, ax = plt.subplots(constrained_layout=True,figsize=(8, 4))
|
||
|
# labels2_h1 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h1 = [PEVDischGrid[0],PEVDischHome[0], PVGrid[0],PVHome[0],PVBattery[0], PESSDischGrid[0],PESSDischHome[0]]
|
||
|
# bottom2_h1 = np.zeros(len(time))
|
||
|
# for label2_h1, data2_h1 in zip(labels2_h1, data_list2_h1):
|
||
|
# ax.bar(time, data2_h1, bottom=bottom2_h1, label=label2_h1)
|
||
|
# bottom2_h1 += data2_h1
|
||
|
# # Set plot labels and title
|
||
|
# ax.set_ylabel('Power Export (kW)', fontname='Times New Roman',fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)',fontname='Times New Roman',fontsize=12)
|
||
|
# # ax.set_title('Power Export (Prosumer-1)', fontname='Times New Roman',fontsize=12)
|
||
|
# ax.legend(loc='upper left',ncol=3,fontsize=12)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 96, 4))
|
||
|
# # Show plot
|
||
|
# plt.show()
|
||
|
|
||
|
# fig_one, ax = plt.subplots(constrained_layout=True, figsize=(9, 6))
|
||
|
#
|
||
|
# labels_h1 = ('EV_Charging', 'Power_EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list_h1 = [PEVCharge[0], PEWH[0], baseload[0], PESSChargeGrid[0]]
|
||
|
# bottom_h1 = np.zeros(len(time))
|
||
|
# colors_h1 = ['red', 'blue', 'green', 'purple']
|
||
|
# for label_h1, data_h1, color_h1 in zip(labels_h1, data_list_h1, colors_h1):
|
||
|
# ax.bar(time, data_h1, bottom=bottom_h1, label=label_h1, color=color_h1)
|
||
|
# bottom_h1 += data_h1
|
||
|
#
|
||
|
# labels2_h1 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h1 = [PEVDischGrid[0], PEVDischHome[0], PVGrid[0], PVHome[0], PVBattery[0], PESSDischGrid[0], PESSDischHome[0]]
|
||
|
# bottom2_h1 = np.zeros(len(time))
|
||
|
# colors2_h1 = ['orange', 'pink', 'brown', 'gray', 'cyan', 'magenta', 'skyblue']
|
||
|
# for label2_h1, data2_h1, color2_h1 in zip(labels2_h1, data_list2_h1, colors2_h1):
|
||
|
# ax.bar(time, np.negative(data2_h1), bottom=bottom2_h1[-len(data2_h1):], label=label2_h1, color=color2_h1) # invert the y-axis by negating the data and bottom values
|
||
|
# bottom2_h1[-len(data2_h1):] -= data2_h1
|
||
|
#
|
||
|
# ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.legend(loc='upper left', fontsize=12, ncol=2)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
#
|
||
|
# plt.show()
|
||
|
#
|
||
|
# my_palette = ['#FF5733', '#C70039', '#900C3F', '#581845', '#2E86C1', '#28B463', '#F7DC6F', '#F39C12', '#D35400', '#8E44AD', '#5D6D7E', '#34495E', '#2C3E50', '#ECF0F1']
|
||
|
#
|
||
|
# palette = sns.color_palette(my_palette)
|
||
|
# fig_one, ax = plt.subplots(constrained_layout=True, figsize=(9, 6))
|
||
|
# labels_h1 = ('EV_Charging', 'Power_EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list_h1 = [PEVCharge[0], PEWH[0], baseload[0], PESSChargeGrid[0]]
|
||
|
# bottom_h1 = np.zeros(len(time))
|
||
|
# for label_h1, data_h1, color_h1 in zip(labels_h1, data_list_h1, palette[:4]):
|
||
|
# ax.bar(time, data_h1, bottom=bottom_h1, label=label_h1, color=color_h1)
|
||
|
# bottom_h1 += data_h1
|
||
|
# labels2_h1 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h1 = [PEVDischGrid[0], PEVDischHome[0], PVGrid[0], PVHome[0], PVBattery[0], PESSDischGrid[0], PESSDischHome[0]]
|
||
|
# bottom2_h1 = np.zeros(len(time))
|
||
|
# for label2_h1, data2_h1, color2_h1 in zip(labels2_h1, data_list2_h1, palette[4:]):
|
||
|
# ax.bar(time, np.negative(data2_h1), bottom=bottom2_h1[-len(data2_h1):], label=label2_h1, color=color2_h1) # invert the y-axis by negating the data and bottom values
|
||
|
# bottom2_h1[-len(data2_h1):] -= data2_h1
|
||
|
# ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.legend(loc='upper left', fontsize=12, ncol=2)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
# def pos_int(x, pos):
|
||
|
# return abs(int(x))
|
||
|
#
|
||
|
# # Format y-axis ticks and labels
|
||
|
# y_ticks = ticker.FuncFormatter(pos_int)
|
||
|
# ax.yaxis.set_major_formatter(y_ticks)
|
||
|
# plt.show()
|
||
|
|
||
|
|
||
|
# fig_two, ax = plt.subplots(constrained_layout=True, figsize=(9, 6))
|
||
|
# labels_h1 = ('EV_Charging', 'Power_EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list_h1 = [PEVCharge[1], PEWH[1], baseload[1], PESSChargeGrid[1]]
|
||
|
# bottom_h1 = np.zeros(len(time))
|
||
|
# for label_h1, data_h1, color_h1 in zip(labels_h1, data_list_h1, palette[:4]):
|
||
|
# ax.bar(time, data_h1, bottom=bottom_h1, label=label_h1, color=color_h1)
|
||
|
# bottom_h1 += data_h1
|
||
|
# labels2_h1 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h1 = [PEVDischGrid[1], PEVDischHome[1], PVGrid[1], PVHome[1], PVBattery[1], PESSDischGrid[1], PESSDischHome[1]]
|
||
|
# bottom2_h1 = np.zeros(len(time))
|
||
|
# for label2_h1, data2_h1, color2_h1 in zip(labels2_h1, data_list2_h1, palette[4:]):
|
||
|
# ax.bar(time, np.negative(data2_h1), bottom=bottom2_h1[-len(data2_h1):], label=label2_h1, color=color2_h1) # invert the y-axis by negating the data and bottom values
|
||
|
# bottom2_h1[-len(data2_h1):] -= data2_h1
|
||
|
# ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.legend(loc='upper left', fontsize=12, ncol=2)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
# def pos_int(x, pos):
|
||
|
# return abs(int(x))
|
||
|
#
|
||
|
# # Format y-axis ticks and labels
|
||
|
# y_ticks = ticker.FuncFormatter(pos_int)
|
||
|
# ax.yaxis.set_major_formatter(y_ticks)
|
||
|
# plt.show()
|
||
|
|
||
|
|
||
|
# fig_three, ax = plt.subplots(constrained_layout=True, figsize=(9, 6))
|
||
|
# labels_h1 = ('EV_Charging', 'Power_EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list_h1 = [PEVCharge[2], PEWH[2], baseload[2], PESSChargeGrid[2]]
|
||
|
# bottom_h1 = np.zeros(len(time))
|
||
|
# for label_h1, data_h1, color_h1 in zip(labels_h1, data_list_h1, palette[:4]):
|
||
|
# ax.bar(time, data_h1, bottom=bottom_h1, label=label_h1, color=color_h1)
|
||
|
# bottom_h1 += data_h1
|
||
|
# labels2_h1 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h1 = [PEVDischGrid[2], PEVDischHome[2], PVGrid[2], PVHome[2], PVBattery[2], PESSDischGrid[2], PESSDischHome[2]]
|
||
|
# bottom2_h1 = np.zeros(len(time))
|
||
|
# for label2_h1, data2_h1, color2_h1 in zip(labels2_h1, data_list2_h1, palette[4:]):
|
||
|
# ax.bar(time, np.negative(data2_h1), bottom=bottom2_h1[-len(data2_h1):], label=label2_h1, color=color2_h1) # invert the y-axis by negating the data and bottom values
|
||
|
# bottom2_h1[-len(data2_h1):] -= data2_h1
|
||
|
# ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.legend(loc='upper left', fontsize=12, ncol=2)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
# def pos_int(x, pos):
|
||
|
# return abs(int(x))
|
||
|
#
|
||
|
# # Format y-axis ticks and labels
|
||
|
# y_ticks = ticker.FuncFormatter(pos_int)
|
||
|
# ax.yaxis.set_major_formatter(y_ticks)
|
||
|
# plt.show()
|
||
|
|
||
|
# fig_four, ax = plt.subplots(constrained_layout=True, figsize=(9, 6))
|
||
|
# labels_h1 = ('EV_Charging', 'Power_EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list_h1 = [PEVCharge[3], PEWH[3], baseload[3], PESSChargeGrid[3]]
|
||
|
# bottom_h1 = np.zeros(len(time))
|
||
|
# for label_h1, data_h1, color_h1 in zip(labels_h1, data_list_h1, palette[:4]):
|
||
|
# ax.bar(time, data_h1, bottom=bottom_h1, label=label_h1, color=color_h1)
|
||
|
# bottom_h1 += data_h1
|
||
|
# labels2_h1 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h1 = [PEVDischGrid[3], PEVDischHome[3], PVGrid[3], PVHome[3], PVBattery[3], PESSDischGrid[3], PESSDischHome[3]]
|
||
|
# bottom2_h1 = np.zeros(len(time))
|
||
|
# for label2_h1, data2_h1, color2_h1 in zip(labels2_h1, data_list2_h1, palette[4:]):
|
||
|
# ax.bar(time, np.negative(data2_h1), bottom=bottom2_h1[-len(data2_h1):], label=label2_h1, color=color2_h1) # invert the y-axis by negating the data and bottom values
|
||
|
# bottom2_h1[-len(data2_h1):] -= data2_h1
|
||
|
# ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.legend(loc='upper left', fontsize=12, ncol=2)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
# def pos_int(x, pos):
|
||
|
# return abs(int(x))
|
||
|
#
|
||
|
# # Format y-axis ticks and labels
|
||
|
# y_ticks = ticker.FuncFormatter(pos_int)
|
||
|
# ax.yaxis.set_major_formatter(y_ticks)
|
||
|
# plt.show()
|
||
|
|
||
|
|
||
|
#
|
||
|
# fig_one, ax = plt.subplots(constrained_layout=True, figsize=(9, 6))
|
||
|
# labels_h1 = ('EV_Charging', 'Power_EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list_h1 = [PEVCharge[0], PEWH[0], baseload[0], PESSChargeGrid[0]]
|
||
|
# bottom_h1 = np.zeros(len(time))
|
||
|
# colors_h1 = ['#FF5733', '#2471A3', '#4CAF50', '#9B59B6']
|
||
|
# for label_h1, data_h1, color_h1 in zip(labels_h1, data_list_h1, colors_h1):
|
||
|
# ax.bar(time, data_h1, bottom=bottom_h1, label=label_h1, color=color_h1)
|
||
|
# bottom_h1 += data_h1
|
||
|
# labels2_h1 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h1 = [PEVDischGrid[0], PEVDischHome[0], PVGrid[0], PVHome[0], PVBattery[0], PESSDischGrid[0], PESSDischHome[0]]
|
||
|
# bottom2_h1 = np.zeros(len(time))
|
||
|
# colors2_h1 = ['#FF8C42', '#3498DB', '#2ECC71', '#8E44AD', '#FFC300', '#E74C3C', '#85C1E9']
|
||
|
# for label2_h1, data2_h1, color2_h1 in zip(labels2_h1, data_list2_h1, colors2_h1):
|
||
|
# ax.bar(time, np.negative(data2_h1), bottom=bottom2_h1[-len(data2_h1):], label=label2_h1, color=color2_h1)
|
||
|
# bottom2_h1[-len(data2_h1):] -= data2_h1
|
||
|
# ax.set_ylabel('Power (kW)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)', fontname='Times New Roman', fontsize=12)
|
||
|
# ax.legend(loc='upper left', fontsize=12, ncol=2)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 98, 4))
|
||
|
# plt.show()
|
||
|
#
|
||
|
|
||
|
#
|
||
|
#
|
||
|
#
|
||
|
#
|
||
|
# fig2_h2, ax = plt.subplots(constrained_layout=True,figsize=(8, 4))
|
||
|
# labels2_h2 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h2 = [PEVDischGrid[1],PEVDischHome[1], PVGrid[1],PVHome[1],PVBattery[1], PESSDischGrid[1],PESSDischHome[1]]
|
||
|
# bottom2_h2 = np.zeros(len(time))
|
||
|
# for label2_h2, data2_h2 in zip(labels2_h2, data_list2_h2):
|
||
|
# ax.bar(time, data2_h2, bottom=bottom2_h2, label=label2_h2)
|
||
|
# bottom2_h2 += data2_h2
|
||
|
# # Set plot labels and title
|
||
|
# ax.set_ylabel('Power Export (kW)', fontname='Times New Roman',fontsize=12)
|
||
|
# ax.set_xlabel('Time (15 min interval)',fontname='Times New Roman',fontsize=12)
|
||
|
# # ax.set_title('Power Export (Prosumer-2)', fontname='Times New Roman',fontsize=12)
|
||
|
# ax.legend(loc='upper left',ncol=3,fontsize=12)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 96, 4))
|
||
|
# # Show plot
|
||
|
# plt.show()
|
||
|
#
|
||
|
#
|
||
|
# fig2_h3, ax = plt.subplots(constrained_layout=True,figsize=(8, 4))
|
||
|
# labels2_h3 = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list2_h3 = [PEVDischGrid[2],PEVDischHome[2], PVGrid[2],PVHome[2],PVBattery[2], PESSDischGrid[2],PESSDischHome[2]]
|
||
|
# bottom2_h3 = np.zeros(len(time))
|
||
|
# for label2_h3, data2_h3 in zip(labels2_h3, data_list2_h3):
|
||
|
# ax.bar(time, data2_h3, bottom=bottom2_h3, label=label2_h3)
|
||
|
# bottom2_h3 += data2_h3
|
||
|
# # Set plot labels and title
|
||
|
# ax.set_ylabel('Power Export (kW)', fontname='Times New Roman')
|
||
|
# ax.set_xlabel('Time (15 min interval)',fontsize=13, fontname='Times New Roman')
|
||
|
# # ax.set_title('Power Export (Consumer-3)', fontname='Times New Roman')
|
||
|
# ax.legend(loc='upper left',ncol=3)
|
||
|
# ax.set_xlim([0, 97])
|
||
|
# ax.set_xticks(np.arange(0, 96, 4))
|
||
|
# # Show plot
|
||
|
# plt.show()
|
||
|
|
||
|
|
||
|
# <editor-fold desc="Home Model for Minimizing Energy">
|
||
|
m2 = py.ConcreteModel()
|
||
|
|
||
|
# <editor-fold desc="Set and Parameters of HOME">
|
||
|
# m2.h1 = py.Set(initialize=sheet_home.home_bus)
|
||
|
#
|
||
|
# m2.h2 = py.Set(initialize=sheet_home.home)
|
||
|
# m2.h = py.Set(m2.h1, m2.h2)
|
||
|
#
|
||
|
#
|
||
|
# def homeset(m):
|
||
|
# return
|
||
|
# # ({h1: h2} for h1 in m2.h1 for h2 in m2.h2)
|
||
|
# # create sets with NO indexes
|
||
|
# m2.setHomes = py.Set(initialize=sheet_home.home.unique())
|
||
|
# m2.setBuses = py.Set(initialize=sheet_home.home_bus.unique())
|
||
|
#
|
||
|
# # create sets WITH indexes
|
||
|
# dict_home_givenBus = {j: sheet_home.home[sheet_home.home_bus == j].unique() for j in m2.setBuses}
|
||
|
# m2.setHome_givenBus = py.Set(m2.setBuses, initialize=dict_home_givenBus)
|
||
|
#
|
||
|
# dict_bus_givenHome = {h:int(sheet_home.home_bus[sheet_home.home==h]) for h in m2.setHomes}
|
||
|
# m2.setBuses_givenHome = py.Param(m2.setHomes, initialize=dict_bus_givenHome)
|
||
|
m2.b = py.Set(initialize=sheet_bus.bus, doc='bus')
|
||
|
|
||
|
m2.h = py.Set(initialize=sheet_home.home)
|
||
|
m2.alpha = py.Param(m2.h, initialize=dict(zip(list(itertools.product(m2.h.data())),
|
||
|
[i for x in sheet_EV.columns if "arrival" in x for i
|
||
|
in
|
||
|
sheet_EV.loc[:, x].values])), doc="Arrival time of EV")
|
||
|
|
||
|
m2.beta = py.Param(m2.h, initialize=dict(zip(list(itertools.product(m2.h.data())),
|
||
|
[i for x in sheet_EV.columns if "departure" in x for i
|
||
|
in
|
||
|
sheet_EV.loc[:, x].values])), doc="departure time of EV")
|
||
|
m2.t_final = py.Param(initialize=96, doc='final time of the day ') # 96
|
||
|
m2.second = py.Param(initialize=2, doc='second time of the day ')
|
||
|
m2.t = py.Set(initialize=sheet_data.time, ordered=True, doc='time period')
|
||
|
# m2.EV_Dur = py.Set(initialize=py.RangeSet(m2.alpha.value, m2.beta.value),
|
||
|
# doc='time interval for EV')
|
||
|
# m2.EV_Dur1 = py.Set(initialize=py.RangeSet(m2.alpha.value + 1, m2.beta.value),
|
||
|
# doc='time interval for EV after arrival')
|
||
|
m2.t_second = py.Set(initialize=py.RangeSet(m2.second.value, m2.t_final.value),
|
||
|
doc='time greater then 1')
|
||
|
m2.t_first = py.Param(initialize=1, doc='first time of the day ')
|
||
|
# m2.ist_interval = py.Set(initialize=py.RangeSet(m2.t_first.value, m2.alpha.value - 1),
|
||
|
# doc='time period before the arrival time')
|
||
|
# m2.second_interval = py.Set(initialize=py.RangeSet(m2.beta.value + 1, m2.t_final.value),
|
||
|
# doc='time interval after departure')
|
||
|
m2.EV_Dur1 = py.Set(m2.h, initialize=lambda m, h: range(m.alpha[h] + 1, m.beta[h] + 1), doc='time interval for EV')
|
||
|
m2.EV_Dur = py.Set(m2.h, initialize=lambda m, h: range(m.alpha[h], m.beta[h] + 1), doc='time interval for EV')
|
||
|
m2.ist_interval = py.Set(m2.h, initialize=lambda m, h: range(m2.t_first.value, m.alpha[h]), doc='time interval for EV')
|
||
|
m2.second_interval = py.Set(m2.h, initialize=lambda m, h: range(m.beta[h] + 1, m2.t_final.value + 1),
|
||
|
doc='time interval after departure')
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="base load and other Parameters of Home">
|
||
|
# Base load : uncontrollable loads
|
||
|
m2.base_load = py.Param(m2.h, m2.t,
|
||
|
initialize=dict(zip(list(itertools.product(m2.h.data(), m2.t.data())),
|
||
|
[i for x in sheet_load.columns if "load" in x for i in
|
||
|
sheet_load.loc[:, x].values])), doc="Base Load")
|
||
|
|
||
|
m2.PV = py.Param(m2.h, m2.t, initialize=dict(zip(list(itertools.product(m2.h.data(), m2.t.data())),
|
||
|
[i for x in sheet_PV.columns if "PV" in x for i in
|
||
|
sheet_PV.loc[:, x].values])), doc="PV Production")
|
||
|
m2.max = py.Param(initialize=10000, doc='maximum value for selling and buying power')
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Electric Vehicle Parameters">
|
||
|
|
||
|
m2.ChargeRate_EV = py.Param(m2.h, initialize=dict(zip(list(itertools.product(m2.h.data())),
|
||
|
[i for x in sheet_EV.columns if "charging_rate" in x for i
|
||
|
in
|
||
|
sheet_EV.loc[:, x].values])), doc="charging rate of EV")
|
||
|
|
||
|
m2.DischRate_EV = py.Param(m2.h, initialize=dict(zip(list(itertools.product(m2.h.data())),
|
||
|
[i for x in sheet_EV.columns if "rate_of_discharging" in x
|
||
|
for i in
|
||
|
sheet_EV.loc[:, x].values])), doc="Discharging rate of EV")
|
||
|
|
||
|
m2.Ch_Effc_EV = py.Param(m2.h, initialize=dict(zip(list(itertools.product(m2.h.data())),
|
||
|
[i for x in sheet_EV.columns if "charging_efficiency" in x for
|
||
|
i in
|
||
|
sheet_EV.loc[:, x].values])), doc="charging efficency of EV")
|
||
|
|
||
|
m2.DischEffc_EV = py.Param(m2.h, initialize=dict(zip(list(itertools.product(m2.h.data())),
|
||
|
[i for x in sheet_EV.columns if
|
||
|
"efficiency_of_discharging" in x for i
|
||
|
in sheet_EV.loc[:, x].values])),
|
||
|
doc="Discharging efficency of EV")
|
||
|
|
||
|
m2.Cap_EV = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EV.columns if "capacity" in x for i in
|
||
|
sheet_EV.loc[:, x].values])), doc="Capacity of EV")
|
||
|
m2.End_percentage_EV = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EV.columns if "end" in x
|
||
|
for i in
|
||
|
sheet_EV.loc[:, x].values])),
|
||
|
doc="Departure energy of EV")
|
||
|
|
||
|
m2.In_Percentage_EV = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EV.columns if "initial" in x
|
||
|
for i in
|
||
|
sheet_EV.loc[:, x].values])),
|
||
|
doc="Initial energy of EV")
|
||
|
|
||
|
m2.Energy_EV_dep = py.Param(m2.h, initialize=dict(
|
||
|
zip(m2.h.data(), np.array(m2.Cap_EV.values()) * m2.End_percentage_EV.values())), # just changing np>py
|
||
|
doc="Departure temperature of EV")
|
||
|
|
||
|
m2.Energy_EV_In = py.Param(m2.h, initialize=dict(
|
||
|
zip(m2.h.data(), np.array(m2.Cap_EV.values()) * m2.In_Percentage_EV.values())),
|
||
|
doc="Initial energy of EV")
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Energy Storage System Parameter">
|
||
|
|
||
|
m2.ChargeRate_ESS = py.Param(m2.h, initialize=dict(zip(list(itertools.product(m2.h.data())),
|
||
|
[i for x in sheet_ESS.columns if "charging_rate" in x for i
|
||
|
in
|
||
|
sheet_ESS.loc[:, x].values])), doc="charging rate of ESS")
|
||
|
|
||
|
m2.DischRate_ESS = py.Param(m2.h, initialize=dict(zip(list(itertools.product(m2.h.data())),
|
||
|
[i for x in sheet_ESS.columns if "rate_of_discharging" in x
|
||
|
for
|
||
|
i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="Discharging rate of ESS")
|
||
|
|
||
|
# m2.ChargeRate_ESS = py.Param(initialize=float(sheet_ESS.charging_rate1), doc='Charging rate of ESS ')
|
||
|
# m2.DischRate_ESS = py.Param(initialize=float(sheet_ESS.discharging_rate1),
|
||
|
# doc='Discharging rate of ESS ')
|
||
|
m2.Cap_ESS = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_ESS.columns if "capacity" in x for i in
|
||
|
sheet_ESS.loc[:, x].values])), doc="Capacity of ESS")
|
||
|
|
||
|
m2.End_percentage_ESS = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_ESS.columns if "end" in x
|
||
|
for i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="Departure energy of ESS")
|
||
|
|
||
|
m2.In_Percentage_ESS = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_ESS.columns if
|
||
|
"initial" in x for i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="Initial energy of ESS")
|
||
|
|
||
|
m2.Ch_Effc_ESS = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_ESS.columns if "charging_efficiency" in x
|
||
|
for i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="charging efficiency of ESS")
|
||
|
|
||
|
m2.DischEffc_ESS = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_ESS.columns if
|
||
|
"efficiency_of_dicharging" in x for
|
||
|
i in
|
||
|
sheet_ESS.loc[:, x].values])),
|
||
|
doc="Discharging efficiency of ESS")
|
||
|
|
||
|
m2.Energy_ESS_In = py.Param(m2.h, initialize=dict(
|
||
|
zip(m2.h.data(), np.array(m2.Cap_ESS.values()) * m2.In_Percentage_ESS.values())),
|
||
|
doc="Initial energy of ESS")
|
||
|
|
||
|
m2.End_En_ESS = py.Param(m2.h, initialize=dict(
|
||
|
zip(m2.h.data(), np.array(m2.Cap_ESS.values()) * m2.End_percentage_ESS.values())),
|
||
|
doc="Departure energy of ESS")
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Electric Water Heater Parameter">
|
||
|
|
||
|
m2.tetta_low = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "tetta_low" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])),
|
||
|
doc="Lower bound of the water temperature")
|
||
|
m2.tetta_up = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "tetta_up" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])),
|
||
|
doc="Upper bound of the water temperature")
|
||
|
m2.tetta_amb_int = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "tetta_amb_init" in x for i
|
||
|
in
|
||
|
sheet_EWH.loc[:, x].values])),
|
||
|
doc="Initial ambient temperature")
|
||
|
|
||
|
m2.tetta_amb = py.Param(m2.h, m2.t,
|
||
|
initialize=dict(zip(list(itertools.product(m2.h.data(), m2.t.data())),
|
||
|
[i for x in sheet_ambient.columns if "celsius" in x for i in
|
||
|
sheet_ambient.loc[:, x].values])), doc="Outdoor temperature")
|
||
|
|
||
|
m2.Q = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "capacity" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Power of the EWH")
|
||
|
|
||
|
m2.R = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "thermal_resistance" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Thermal resistance")
|
||
|
m2.C = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "thermal_capacitance" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Thermal resistance")
|
||
|
|
||
|
m2.M = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "water_cap" in x for i in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Water Capacity (L)")
|
||
|
|
||
|
m2.tetta_EWH_int = py.Param(m2.h, initialize=dict(zip((m2.h.data()),
|
||
|
[i for x in sheet_EWH.columns if "tetta_wat_init" in x for i
|
||
|
in
|
||
|
sheet_EWH.loc[:, x].values])), doc="Initial temperature")
|
||
|
|
||
|
m2.water_use = py.Param(m2.h, m2.t,
|
||
|
initialize=dict(zip(list(itertools.product(m2.h.data(), m2.t.data())),
|
||
|
[i for x in sheet_wateruse.columns if "Litre" in x for i in
|
||
|
sheet_wateruse.loc[:, x].values])), doc="Hot Water usage")
|
||
|
|
||
|
# </editor-fold2
|
||
|
|
||
|
# <editor-fold desc="Price, Base load, and time duration parameters">
|
||
|
# m2.Buy_price = py.Param(m2.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_RTP)), doc='Buying Price')
|
||
|
# m2.Buy_price = py.Param(m2.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_TOU_winter)), doc='Buying Price')
|
||
|
m2.Buy_price = py.Param(m2.t, initialize=dict(zip(sheet_data.time, sheet_data.TOU_PGE)), doc='Buying Price')
|
||
|
# m2.Buy_price = py.Param(m2.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_flex_D)), doc='Buying Price')
|
||
|
# m2.Buy_price = py.Param(m2.t, initialize=dict(zip(sheet_data.time, sheet_data.Buy_price_fixed)), doc='Buying Price')
|
||
|
m2.Sell_price = py.Param(m2.t, initialize=dict(zip(sheet_data.time, sheet_data.FiT)), doc='Selling Price')
|
||
|
# Time duration: we took 15 mint granularity so in one hour it will be 1/4
|
||
|
m2.time_d = py.Param(initialize=(1 / 4), doc='time duration ')
|
||
|
# m2.DPT = py.Param(initialize=0.069, doc='Daily power tariff ') # .069
|
||
|
# </editor-fold>
|
||
|
|
||
|
print('Code Starts for HEMS Energy Optimization')
|
||
|
# Variable
|
||
|
|
||
|
# <editor-fold desc="buying and selling variables">
|
||
|
m2.P_Buy_Grid = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.P_Sell_Grid = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.S_P_sell = py.Var(m2.h, m2.t, within=py.Binary)
|
||
|
m2.S_P_buy = py.Var(m2.h, m2.t, within=py.Binary)
|
||
|
m2.peak = py.Var(within=py.NonNegativeReals)
|
||
|
# m2.pflex = py.Var(m2.b, m2.t, doc='active Flexibility')
|
||
|
# m2.qflex = py.Var(m2.b, m2.t,
|
||
|
# doc='reactive Flexibility') # i comment this because i use m2.pflex variable * 0.48
|
||
|
m2.pmax = py.Var(m2.h, m2.t, bounds=(0, None), doc='Maximum power')
|
||
|
m2.pmax_int = py.Var(m2.h, m2.t, bounds=(0, None), doc='Initial Power')
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Electric vehicle Variables">
|
||
|
m2.P_EV_Charge = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.P_EV_Disch = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.P_EV_Disch_Home = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.P_EV_Disch_Grid = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.S_EV_Charge = py.Var(m2.h, m2.t, within=py.Binary)
|
||
|
m2.S_EV_Disch = py.Var(m2.h, m2.t, within=py.Binary)
|
||
|
m2.Energy_EV = py.Var(m2.h, m2.t, bounds=(0, None)) # SOC of the EV
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Energy storage system Variables">
|
||
|
m2.Energy_ESS = py.Var(m2.h, m2.t, bounds=(0, None)) # SOC of the ESS
|
||
|
m2.S_ESS_Charge = py.Var(m2.h, m2.t, within=py.Binary)
|
||
|
m2.P_ESS_Disch = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.P_ESS_Charge = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.P_ESS_Charge_Grid = py.Var(m2.h, m2.t, bounds=(0, None), doc='ESS charging from the Grid')
|
||
|
m2.P_ESS_Disch_Home = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.P_ESS_Disch_Grid = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.S_ESS_Disch = py.Var(m2.h, m2.t, within=py.Binary)
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="PV variables">
|
||
|
m2.PV_Home = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.PV_Grid = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
m2.PV_Battery = py.Var(m2.h, m2.t, bounds=(0, None))
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="EWH Variables">
|
||
|
m2.tetta_EWH_wat = py.Var(m2.h, m2.t)
|
||
|
m2.S_EWH = py.Var(m2.h, m2.t, within=py.Binary)
|
||
|
m2.P_EWH = py.Var(m2.h, m2.t, doc='power of EWH')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
# Constraints
|
||
|
|
||
|
|
||
|
# <editor-fold desc="Buying and Selling Power">
|
||
|
def Power_buy(m2, h, i):
|
||
|
return m2.P_Buy_Grid[h, i] == m2.base_load[h, i] + m2.Q[h] * m2.S_EWH[h, i] + m2.P_EV_Charge[h, i] - \
|
||
|
m2.P_EV_Disch_Home[h, i] + m2.P_ESS_Charge[h, i] - \
|
||
|
m2.P_ESS_Disch_Home[h, i] - m2.PV_Home[h, i] - m2.PV_Battery[h, i]
|
||
|
|
||
|
|
||
|
m2.Const_1 = py.Constraint(m2.h, m2.t, rule=Power_buy, doc='Power buy from the Grid')
|
||
|
|
||
|
|
||
|
def Power_buy2(m2, h, i):
|
||
|
return m2.P_Buy_Grid[h, i] <= m2.max * m2.S_P_buy[h, i]
|
||
|
|
||
|
|
||
|
m2.Const_1a = py.Constraint(m2.h, m2.t, rule=Power_buy2,
|
||
|
doc='removing the nonlinearity in the objective fucntion')
|
||
|
|
||
|
|
||
|
def Power_sell1(m2, h, i):
|
||
|
return m2.P_Sell_Grid[h, i] == m2.P_EV_Disch_Grid[h, i] + m2.PV_Grid[h, i] + m2.P_ESS_Disch_Grid[h, i]
|
||
|
|
||
|
|
||
|
m2.Const_2 = py.Constraint(m2.h, m2.t, rule=Power_sell1, doc='Power sell to the Grid')
|
||
|
|
||
|
|
||
|
# def Power_sell2(m2, h, i):
|
||
|
# return m2.P_Sell_Grid[h, i] <= m2.max * m2.S_P_sell[h, i]
|
||
|
#
|
||
|
#
|
||
|
# m2.Const_2a = py.Constraint(m2.h, m2.t, rule=Power_sell2, doc='Power sell to the Grid')
|
||
|
def Power_sell2(m2, h, i):
|
||
|
return m2.P_Sell_Grid[h, i] <= m2.max.value * (1 - m2.S_P_buy[h, i])
|
||
|
|
||
|
|
||
|
m2.Const_2a = py.Constraint(m2.h, m2.t, rule=Power_sell2, doc='Power sell to the Grid')
|
||
|
|
||
|
|
||
|
# def Status_Power(m2, h, i):
|
||
|
# return m2.S_P_buy[h, i] + m2.S_P_sell[h, i] <= 1
|
||
|
#
|
||
|
#
|
||
|
# m2.Const_3 = py.Constraint(m2.h, m2.t, rule=Status_Power,
|
||
|
# doc='Buying and selling power will not occur at same time')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
|
||
|
# <editor-fold desc="Electric vehicle Constraints">
|
||
|
def Power_EV_Charge_limit1(m2, h, i):
|
||
|
return m2.P_EV_Charge[h, i] <= m2.ChargeRate_EV[h] * m2.S_EV_Charge[h, i]
|
||
|
|
||
|
|
||
|
m2.Const_EV_1 = py.Constraint(m2.h, m2.t, rule=Power_EV_Charge_limit1,
|
||
|
doc='Charging power of EV (Upper limit)')
|
||
|
|
||
|
|
||
|
def Power_EV_Charge_limit2(m2, h, i):
|
||
|
return m2.P_EV_Charge[h, i] >= 0
|
||
|
|
||
|
|
||
|
m2.Const_EV_2 = py.Constraint(m2.h, m2.t, rule=Power_EV_Charge_limit2,
|
||
|
doc='Charging power of EV (lower limit)')
|
||
|
|
||
|
|
||
|
def Power_EV_Disch_limit1(m2, h, i):
|
||
|
return m2.P_EV_Disch[h, i] <= (m2.DischRate_EV[h] * m2.S_EV_Disch[h, i])
|
||
|
|
||
|
|
||
|
m2.Const_EV_3 = py.Constraint(m2.h, m2.t, rule=Power_EV_Disch_limit1,
|
||
|
doc='Discharging power of EV (Upper limit)')
|
||
|
|
||
|
|
||
|
def Power_EV_Disch_limit2(m2, h, i):
|
||
|
return m2.P_EV_Disch[h, i] >= 0
|
||
|
|
||
|
|
||
|
m2.Const_EV_4 = py.Constraint(m2.h, m2.t, rule=Power_EV_Disch_limit2,
|
||
|
doc='Discharging power of EV (Lower limit)')
|
||
|
|
||
|
|
||
|
def Status_EV(m2, h, i):
|
||
|
return m2.S_EV_Disch[h, i] + m2.S_EV_Charge[h, i] <= 1
|
||
|
|
||
|
|
||
|
m2.Const_EV_5 = py.Constraint(m2.h, m2.t, rule=Status_EV,
|
||
|
doc='Charging and discharging will not occur at same time')
|
||
|
|
||
|
|
||
|
def Power_EV_Disch(m2, h, i):
|
||
|
return m2.P_EV_Disch[h, i] * m2.DischEffc_EV[h] == (
|
||
|
(m2.P_EV_Disch_Home[h, i] + m2.P_EV_Disch_Grid[h, i]))
|
||
|
|
||
|
|
||
|
m2.Const_EV_6 = py.Constraint(m2.h, m2.t, rule=Power_EV_Disch, doc='Discharging power of EV to '
|
||
|
'Home and Grid')
|
||
|
|
||
|
|
||
|
def SoC_EV1(m2, h, i):
|
||
|
return m2.Energy_EV[h, i] == m2.Energy_EV_In[h]
|
||
|
|
||
|
|
||
|
# m2.Const_EV_7a = py.Constraint(m2.h, [m2.alpha.value], rule=SoC_EV1, doc='SoC of the EV at arrival')
|
||
|
m2.Const_EV_7a = py.Constraint([(h, py.value(m2.alpha[h])) for h in m2.h], rule=SoC_EV1, doc='SoC of the EV at arrival')
|
||
|
|
||
|
|
||
|
def SoC_EV2(m2, h, i):
|
||
|
return m2.Energy_EV[h, i] == m2.Energy_EV[h, i - 1] + (
|
||
|
m2.P_EV_Charge[h, i] * m2.Ch_Effc_EV[h] - m2.P_EV_Disch[h, i]) * m2.time_d
|
||
|
|
||
|
|
||
|
# m2.Const_EV_7b = py.Constraint(m2.h, m2.EV_Dur1, rule=SoC_EV2, doc='SoC of the EV after arrival time')
|
||
|
m2.Const_EV_7b = py.ConstraintList()
|
||
|
for h in m2.h:
|
||
|
for i in m2.EV_Dur1[h].value:
|
||
|
m2.Const_EV_7b.add(SoC_EV2(m2, h, i))
|
||
|
|
||
|
|
||
|
def EV_availability1(m2, h, i):
|
||
|
return m2.Energy_EV[h, i] == 0
|
||
|
|
||
|
|
||
|
# m2.Const_EV_8a = py.Constraint(m2.h, m2.ist_interval, rule=EV_availability1,
|
||
|
# doc='SOC available before arrival time')
|
||
|
m2.Const_EV_8a = py.ConstraintList()
|
||
|
for h in m2.h:
|
||
|
for i in m2.ist_interval[h].value:
|
||
|
m2.Const_EV_8a.add(EV_availability1(m2, h, i))
|
||
|
|
||
|
|
||
|
def EV_availability2(m2, h, i):
|
||
|
return m2.Energy_EV[h, i] == 0
|
||
|
|
||
|
|
||
|
# m2.Const_EV_8b = py.Constraint(m2.h, m2.second_interval, rule=EV_availability2,
|
||
|
# doc='SOC available after departure time')
|
||
|
m2.Const_EV_8b = py.ConstraintList()
|
||
|
for h in m2.h:
|
||
|
for i in m2.second_interval[h].value:
|
||
|
m2.Const_EV_8b.add(EV_availability2(m2, h, i))
|
||
|
|
||
|
|
||
|
def EV_status_available1(m2, h, i):
|
||
|
return m2.S_EV_Disch[h, i] + m2.S_EV_Charge[h, i] == 0
|
||
|
|
||
|
|
||
|
# m2.Const_EV_9a = py.Constraint(m2.h, m2.ist_interval, rule=EV_status_available1,
|
||
|
# doc='EV availability before arrival time')
|
||
|
m2.Const_EV_9a = py.ConstraintList()
|
||
|
for h in m2.h:
|
||
|
for i in m2.ist_interval[h].value:
|
||
|
m2.Const_EV_9a.add(EV_status_available1(m2, h, i))
|
||
|
|
||
|
|
||
|
def EV_status_available2(m2, h, i):
|
||
|
return m2.S_EV_Disch[h, i] + m2.S_EV_Charge[h, i] == 0
|
||
|
|
||
|
|
||
|
# m2.Const_EV_9b = py.Constraint(m2.h, m2.second_interval, rule=EV_status_available2,
|
||
|
# doc='EV availability after arrival time')
|
||
|
m2.Const_EV_9b = py.ConstraintList()
|
||
|
for h in m2.h:
|
||
|
for i in m2.second_interval[h].value:
|
||
|
m2.Const_EV_9b.add(EV_status_available2(m2, h, i))
|
||
|
|
||
|
|
||
|
def EV_SoC_limit1(m2, h, i):
|
||
|
return m2.Energy_EV[h, i] >= 0.2 * m2.Cap_EV[h]
|
||
|
|
||
|
|
||
|
# m2.Const_EV_10 = py.Constraint(m2.h, m2.EV_Dur, rule=EV_SoC_limit1, doc='Minimum SoC of EV')
|
||
|
m2.Const_EV_10 = py.ConstraintList()
|
||
|
for h in m2.h:
|
||
|
for i in m2.EV_Dur[h].value:
|
||
|
m2.Const_EV_10.add(EV_SoC_limit1(m2, h, i))
|
||
|
|
||
|
|
||
|
def EV_SoC_limit2(m2, h, i):
|
||
|
return m2.Energy_EV[h, i] <= m2.Cap_EV[h]
|
||
|
|
||
|
|
||
|
# m2.Const_EV_11 = py.Constraint(m2.h, m2.t, rule=EV_SoC_limit2, doc='Maximum SoC of EV')
|
||
|
m2.Const_EV_11 = py.ConstraintList()
|
||
|
for h in m2.h:
|
||
|
for i in m2.EV_Dur[h].value:
|
||
|
m2.Const_EV_11.add(EV_SoC_limit2(m2, h, i))
|
||
|
|
||
|
|
||
|
def EV_final_SoC(m2, h, i):
|
||
|
return m2.Energy_EV[h, i] == m2.Energy_EV_dep[h]
|
||
|
|
||
|
|
||
|
# m2.Const_EV_12 = py.Constraint(m2.h, [m2.beta.value], rule=EV_final_SoC,
|
||
|
# doc='Final SoC of EV at departure time')
|
||
|
m2.Const_EV_12 = py.Constraint([(h, py.value(m2.beta[h])) for h in m2.h], rule=EV_final_SoC,
|
||
|
doc='Final SoC of EV at departure time')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Energy Storage System Constraints">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
||
|
def SoC_ESS1(m2, h, i):
|
||
|
return m2.Energy_ESS[h, i] == m2.Energy_ESS_In[h]
|
||
|
|
||
|
|
||
|
m2.Const_ESS_1a = py.Constraint(m2.h, [m2.t_first.value], rule=SoC_ESS1,
|
||
|
doc='SoC of ESS') # need to check wather we need to put in the brackets or not
|
||
|
|
||
|
|
||
|
def SoC_ESS2(m2, h, i):
|
||
|
return m2.Energy_ESS[h, i] == m2.Energy_ESS[h, i - 1] + (
|
||
|
m2.P_ESS_Charge[h, i] * m2.Ch_Effc_ESS[h] - m2.P_ESS_Disch[h, i]) * m2.time_d
|
||
|
|
||
|
|
||
|
m2.Const_ESS_1b = py.Constraint(m2.h, m2.t_second, rule=SoC_ESS2, doc='SoC of ESS')
|
||
|
|
||
|
|
||
|
def Power_ESS_Charge(m2, h, i):
|
||
|
return m2.P_ESS_Charge[h, i] == m2.P_ESS_Charge_Grid[h, i] + m2.PV_Battery[h, i]
|
||
|
|
||
|
|
||
|
m2.Const_ESS_2 = py.Constraint(m2.h, m2.t, rule=Power_ESS_Charge,
|
||
|
doc='Charging power of ESS ')
|
||
|
|
||
|
|
||
|
def Power_ESS_Charge_limit1(m2, h, i):
|
||
|
return (m2.P_ESS_Charge_Grid[h, i] + m2.PV_Battery[h, i]) <= m2.ChargeRate_ESS[h] * m2.S_ESS_Charge[
|
||
|
h, i]
|
||
|
|
||
|
|
||
|
m2.Const_ESS_2a = py.Constraint(m2.h, m2.t, rule=Power_ESS_Charge_limit1,
|
||
|
doc='Charging power of ESS (Upper limit)')
|
||
|
|
||
|
|
||
|
def Power_ESS_Charge_limit2(m2, h, i):
|
||
|
return (m2.P_ESS_Charge_Grid[h, i] + m2.PV_Battery[h, i]) >= 0
|
||
|
|
||
|
|
||
|
m2.Const_ESS_2b = py.Constraint(m2.h, m2.t, rule=Power_ESS_Charge_limit2,
|
||
|
doc='Charging power of ESS (Lower limit)')
|
||
|
|
||
|
|
||
|
def Power_ESS_Disch_limit1(m2, h, i):
|
||
|
return (m2.P_ESS_Disch_Home[h, i] + m2.P_ESS_Disch_Grid[h, i]) <= (
|
||
|
m2.DischRate_ESS[h] * m2.S_ESS_Disch[h, i])
|
||
|
|
||
|
|
||
|
m2.Const_ESS_4 = py.Constraint(m2.h, m2.t, rule=Power_ESS_Disch_limit1,
|
||
|
doc='Discharging power of ESS (Upper limit)')
|
||
|
|
||
|
|
||
|
def Power_ESS_Disch_limit2(m2, h, i):
|
||
|
return (m2.P_ESS_Disch_Home[h, i] + m2.P_ESS_Disch_Grid[h, i]) >= 0
|
||
|
|
||
|
|
||
|
m2.Const_ESS_5 = py.Constraint(m2.h, m2.t, rule=Power_ESS_Disch_limit2,
|
||
|
doc='Discharging power of ESS (Lower limit)')
|
||
|
|
||
|
|
||
|
def Status_ESS(m2, h, i):
|
||
|
return m2.S_ESS_Disch[h, i] + m2.S_ESS_Charge[h, i] <= 1
|
||
|
|
||
|
|
||
|
m2.Const_ESS_6 = py.Constraint(m2.h, m2.t, rule=Status_ESS,
|
||
|
doc='Charging and discharging will not occur at same time')
|
||
|
|
||
|
|
||
|
def Power_ESS_Disch(m2, h, i):
|
||
|
return (m2.P_ESS_Disch[h, i]) * m2.DischEffc_ESS[h] == (m2.P_ESS_Disch_Home[h, i] + m2.P_ESS_Disch_Grid[h, i])
|
||
|
|
||
|
|
||
|
m2.Const_ESS_7 = py.Constraint(m2.h, m2.t, rule=Power_ESS_Disch,
|
||
|
doc='Discharging power of ESS to Home and Grid')
|
||
|
|
||
|
|
||
|
def ESS_SoC_limit1(m2, h, i):
|
||
|
return m2.Energy_ESS[h, i] >= 0.1 * m2.Cap_ESS[h]
|
||
|
|
||
|
|
||
|
m2.Const_ESS_8 = py.Constraint(m2.h, m2.t, rule=ESS_SoC_limit1, doc='Minimum SoC of ESS')
|
||
|
|
||
|
|
||
|
def ESS_SoC_limit2(m2, h, i):
|
||
|
return m2.Energy_ESS[h, i] <= m2.Cap_ESS[h]
|
||
|
|
||
|
|
||
|
m2.Const_ESS_9 = py.Constraint(m2.h, m2.t, rule=ESS_SoC_limit2, doc='Maximum SoC of ESS')
|
||
|
|
||
|
|
||
|
def ESS_final_SoC(m2, h, i):
|
||
|
return m2.Energy_ESS[h, i] >= m2.End_En_ESS[h]
|
||
|
|
||
|
|
||
|
m2.Const_ESS_10 = py.Constraint(m2.h, [m2.t_final.value], rule=ESS_final_SoC,
|
||
|
doc='Final SoC of ESS at departure time')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="Electric Water Heater (EWH) Constraints">
|
||
|
def EWH_limit1(m2, h, i):
|
||
|
return m2.tetta_EWH_wat[h, i] <= m2.tetta_up[h]
|
||
|
|
||
|
|
||
|
m2.Const_EWH_1 = py.Constraint(m2.h, m2.t, rule=EWH_limit1, doc='Maximum Limit')
|
||
|
|
||
|
|
||
|
def EWH_limit2(m2, h, i):
|
||
|
return m2.tetta_EWH_wat[h, i] >= m2.tetta_low[h]
|
||
|
|
||
|
|
||
|
m2.Const_EWH_2 = py.Constraint(m2.h, m2.t, rule=EWH_limit2, doc='Minimum Limit')
|
||
|
|
||
|
|
||
|
def EWH_power(m2, h, i):
|
||
|
return m2.P_EWH[h, i] == m2.Q[h] * m2.S_EWH[h, i]
|
||
|
|
||
|
|
||
|
m2.Const_EWH_3 = py.Constraint(m2.h, m2.t, rule=EWH_power, doc='Electric water heater power')
|
||
|
|
||
|
|
||
|
def EWH_temp_1(m2, h, i):
|
||
|
return m2.tetta_EWH_wat[h, i] == m2.tetta_amb[h, i] + m2.Q[h] * m2.S_EWH[h, i] * m2.R[
|
||
|
h] * m2.time_d - (
|
||
|
(m2.M[h] - m2.water_use[h, i]) / m2.M[h]) * (
|
||
|
m2.tetta_amb_int[h] - m2.tetta_EWH_int[h]) * py.exp(-m2.time_d / (m2.R[h] * m2.C[h]))
|
||
|
|
||
|
|
||
|
m2.Const_EWH_4 = py.Constraint(m2.h, [m2.t_first.value], rule=EWH_temp_1, doc='EWH m2')
|
||
|
|
||
|
|
||
|
def EWH_temp_2(m2, h, i):
|
||
|
return m2.tetta_EWH_wat[h, i] == m2.tetta_amb[h, i] + m2.Q[h] * m2.S_EWH[h, i] * m2.R[
|
||
|
h] * m2.time_d - (
|
||
|
(m2.M[h] - m2.water_use[h, i]) / m2.M[h]) * (
|
||
|
m2.tetta_amb[h, i] - m2.tetta_EWH_wat[h, i - 1]) * py.exp(
|
||
|
-m2.time_d / (m2.R[h] * m2.C[h]))
|
||
|
|
||
|
|
||
|
m2.Const_EWH_5 = py.Constraint(m2.h, m2.t_second, rule=EWH_temp_2, doc='EWH m2')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
# <editor-fold desc="PV Constraints"> # need to include the PV equation of omer paper
|
||
|
|
||
|
def PV_production(m2, h, i):
|
||
|
return m2.PV_Grid[h, i] + m2.PV_Home[h, i] + m2.PV_Battery[h, i] == m2.PV[h, i]
|
||
|
|
||
|
|
||
|
m2.Const_PV = py.Constraint(m2.h, m2.t, rule=PV_production, doc='PV Production')
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
# <editor-fold desc="objective functions">
|
||
|
|
||
|
|
||
|
def objective_rule(m2):
|
||
|
return sum(m2.P_Buy_Grid[h, i] for i in m2.t for h in m2.h) * m2.time_d.value
|
||
|
|
||
|
|
||
|
m2.obj2 = py.Objective(rule=objective_rule, sense=py.minimize, doc='Definition of objective function')
|
||
|
# <editor-fold desc="Solver">
|
||
|
m2.write('m22.lp', io_options={'symbolic_solver_labels': True})
|
||
|
opt = py.SolverFactory('gurobi')
|
||
|
opt.options["mipgap"] = 0.8 # 0.155 for load and 0.8 for other load
|
||
|
result = opt.solve(m2, tee=True) # report_timing=True
|
||
|
|
||
|
print(result)
|
||
|
# print('Sum of all homes objectives = ', py.value(m2.obj2))
|
||
|
# for h in m2.h:
|
||
|
# print(f"Objective of Home {h} : " f"{(sum(m2.P_Buy_Grid[h, i].value for i in m2.t) * m2.time_d.value)}")
|
||
|
# print("second objective Done ")
|
||
|
|
||
|
|
||
|
# </editor-fold>
|
||
|
# # <editor-fold desc="Ploting">
|
||
|
|
||
|
#
|
||
|
# # # <editor-fold desc="Ploting">
|
||
|
# #
|
||
|
# # # Plotting for automatic ploting
|
||
|
time = [i for i in m2.t]
|
||
|
pl2_pbuy = [[] for j in m2.h]
|
||
|
# pl2_psell = [[] for j in m2.h]
|
||
|
PEVCharge2 = [[] for j in m2.h]
|
||
|
PESSChargeGrid2 = [[] for j in m2.h]
|
||
|
PEVDischHome2 = [[] for j in m2.h]
|
||
|
PESSDischHome2 = [[] for j in m2.h]
|
||
|
# # # SEVCharge = [[] for j in m2.h]
|
||
|
# # # SEVDisch = [[] for j in m2.h]
|
||
|
PEVDischGrid2 = [[] for j in m2.h]
|
||
|
PESSDischGrid2 = [[] for j in m2.h]
|
||
|
PVHome2 = [[] for j in m2.h]
|
||
|
PVGrid2 = [[] for j in m2.h]
|
||
|
PVBattery2 = [[] for j in m2.h]
|
||
|
# # # PPV = [[] for j in m2.h]
|
||
|
baseload2 = [[] for j in m2.h]
|
||
|
# # # EnergyESS = [[] for j in m2.h]
|
||
|
# # # EnergyEV = [[] for j in m2.h]
|
||
|
# # # tettaEWHwat = [[] for j in m2.h]
|
||
|
# # # SEWH = [[] for j in m2.h]
|
||
|
PEWH2 = [[] for j in m2.h]
|
||
|
# # # Buyprice = [] # cost = [i for i in m2.c.values()] <<< check this
|
||
|
# # # Sellprice = []
|
||
|
# # #
|
||
|
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.P_Buy_Grid.items():
|
||
|
if k[0] == j:
|
||
|
pl2_pbuy[j - 1].append(py.value(v))
|
||
|
# for j in m2.h:
|
||
|
# for k, v in m2.P_Sell_Grid.items():
|
||
|
# if k[0] == j:
|
||
|
# pl2_psell[j - 1].append(py.value(v))
|
||
|
# # #
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.P_EV_Charge.items():
|
||
|
if k[0] == j:
|
||
|
PEVCharge2[j - 1].append(py.value(v))
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.P_EV_Disch_Home.items():
|
||
|
if k[0] == j:
|
||
|
PEVDischHome2[j - 1].append(py.value(v))
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.P_EV_Disch_Grid.items():
|
||
|
if k[0] == j:
|
||
|
PEVDischGrid2[j - 1].append(py.value(v))
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.P_ESS_Charge_Grid.items():
|
||
|
if k[0] == j:
|
||
|
PESSChargeGrid2[j - 1].append(py.value(v))
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.P_ESS_Disch_Home.items():
|
||
|
if k[0] == j:
|
||
|
PESSDischHome2[j - 1].append(py.value(v))
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.P_ESS_Disch_Grid.items():
|
||
|
if k[0] == j:
|
||
|
PESSDischGrid2[j - 1].append(py.value(v))
|
||
|
# # #
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.PV_Home.items():
|
||
|
if k[0] == j:
|
||
|
PVHome2[j - 1].append(py.value(v))
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.PV_Grid.items():
|
||
|
if k[0] == j:
|
||
|
PVGrid2[j - 1].append(py.value(v))
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.PV_Battery.items():
|
||
|
if k[0] == j:
|
||
|
PVBattery2[j - 1].append(py.value(v))
|
||
|
# # # for j in m2.h:
|
||
|
# # # for k, v in m2.PV.items():
|
||
|
# # # if k[0] == j:
|
||
|
# # # PPV[j - 1].append(py.value(v))
|
||
|
# # #
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.base_load.items():
|
||
|
if k[0] == j:
|
||
|
baseload2[j - 1].append(py.value(v))
|
||
|
# # #
|
||
|
# # # for j in m2.h:
|
||
|
# # # for k, v in m2.Energy_ESS.items():
|
||
|
# # # if k[0] == j:
|
||
|
# # # EnergyESS[j - 1].append(py.value(v))
|
||
|
# # # for j in m2.h:
|
||
|
# # # for k, v in m2.Energy_EV.items():
|
||
|
# # # if k[0] == j:
|
||
|
# # # EnergyEV[j - 1].append(py.value(v))
|
||
|
for j in m2.h:
|
||
|
for k, v in m2.P_EWH.items():
|
||
|
if k[0] == j:
|
||
|
PEWH2[j - 1].append(py.value(v))
|
||
|
# # # for j in m2.h:
|
||
|
# # # for k, v in m2.tetta_EWH_wat.items():
|
||
|
# # # if k[0] == j:
|
||
|
# # # tettaEWHwat[j - 1].append(py.value(v))
|
||
|
# # #
|
||
|
# # # for j in m2.h:
|
||
|
# # # for k, v in m2.S_EWH.items():
|
||
|
# # # if k[0] == j:
|
||
|
# # # SEWH[j - 1].append(py.value(v))
|
||
|
# # #
|
||
|
|
||
|
# # #
|
||
|
# # # # alone
|
||
|
# # for k in range(len(m2.h)):
|
||
|
# # fig, ax = plt.subplots(5, 2, figsize=(10, 10))
|
||
|
# # ax[0, 0].bar(time, pl2_pbuy[k], label='Buying power-obj2')
|
||
|
# # ax[0, 0].bar(time, pl1_pbuy[k], label='Buying power-obj1')
|
||
|
# # ax[1, 0].bar(time, pl2_psell[k], label='Selling power-obj2')
|
||
|
# # ax[1, 0].bar(time, pl1_psell[k], label='Selling power-obj1')
|
||
|
# # ax[0, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # ax[1, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[0, 1].bar(time, baseload[k], label='Base load', color='r')
|
||
|
# # # # ax[0, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[0, 2].plot(time, Buyprice, label='Buyprice')
|
||
|
# # # # ax[0, 2].plot(time, Sellprice, label='Sellprice')
|
||
|
# # # # ax[0, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[1, 0].bar(time, PEWH[k], label='EWH Power', color='r')
|
||
|
# # # # ax[1, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[1, 1].plot(time, tettaEWHwat[k], label='EWH Temp')
|
||
|
# # # # ax[1, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[1, 2].bar(time, SEWH[k], label='Status of EWH')
|
||
|
# # # # ax[1, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[2, 0].bar(time, PESSChargeGrid[k], label='ESS Charging from Grid', color='r')
|
||
|
# # # # ax[2, 0].bar(time, PVBattery[k], label='ESS Charging from PV', color='g')
|
||
|
# # # # ax[2, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[2, 1].bar(time, PESSDischHome[k], label='ESS Disch to home')
|
||
|
# # # # ax[2, 1].bar(time, PESSDischGrid[k], label='ESS Disch to grid')
|
||
|
# # # # ax[2, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[2, 2].plot(time, EnergyESS[k], label='Energy of ESS', color='g')
|
||
|
# # # # ax[2, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[3, 0].bar(time, PVHome[k], label='PV to Home')
|
||
|
# # # # ax[3, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[3, 1].bar(time, PVGrid[k], label='PV to Grid')
|
||
|
# # # # ax[3, 1].bar(time, PVBattery[k], label='PV to battery')
|
||
|
# # # # ax[3, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[3, 2].bar(time, PVHome[k], label='PV to Home')
|
||
|
# # # # ax[3, 2].bar(time, PVGrid[k], label='PV to Grid')
|
||
|
# # # # ax[3, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[4, 0].bar(time, PEVCharge[k], label='EV Charging power', color='r')
|
||
|
# # # # ax[4, 0].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[4, 1].bar(time, PEVDischHome[k], label='EV Disch to home')
|
||
|
# # # # ax[4, 1].bar(time, PEVDischGrid[k], label='EV Disch to grid')
|
||
|
# # # # ax[4, 1].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[4, 2].plot(time, EnergyEV[k], label='Energy of EV', color='r')
|
||
|
# # # # ax[4, 2].legend(loc='best', fontsize='small', ncol=3)
|
||
|
# # # # ax[4, 0].set_xlabel('Time (step)')
|
||
|
# # # # ax[4, 1].set_xlabel('Time (step)')
|
||
|
# # # # ax[4, 2].set_xlabel('Time (step)')
|
||
|
# # # # ax[0, 2].set_ylabel('price ($/W/h)')
|
||
|
# # # # ax[0, 0].set_ylabel('Power1 (kW)')
|
||
|
# # # # ax[1, 0].set_ylabel('Power (kW)')
|
||
|
# # # # ax[2, 0].set_ylabel('Power (kW)')
|
||
|
# # # # ax[3, 0].set_ylabel('Power (kW)')
|
||
|
# # # # ax[4, 0].set_ylabel('Power (kW)')
|
||
|
# # # # ax[2, 2].set_ylabel('Energy (kWh)')
|
||
|
# # # # ax[4, 2].set_ylabel('Energy (kWh)')
|
||
|
# # plt.suptitle(
|
||
|
# # f" Home {k + 1} : " f" Cost is : {py.value((sum(m2.P_Buy_Grid[k + 1, i] for i in m2.t) * m2.time_d))}")
|
||
|
# # pl.tight_layout()
|
||
|
# # plt.show()
|
||
|
# # #
|
||
|
# # # # code for if you want to plot all the home togathoer
|
||
|
# fig, ax = plt.subplots(4, len(m2.h), figsize=(10, 10), constrained_layout=True, sharex='col', sharey='row')
|
||
|
# gs = fig.add_gridspec(hspace=0, wspace=0)
|
||
|
# for k in range(len(m2.h)):
|
||
|
# ax[0, k].bar(time, pl1_pbuy[k], color='green')
|
||
|
# ax[1, k].bar(time, pl2_pbuy[k], color='blue')
|
||
|
# ax[2, k].bar(time, pl1_psell[k], color='green')
|
||
|
# ax[3, k].bar(time, pl2_psell[k], label='blue')
|
||
|
# # ax[0, k].legend(loc='best')
|
||
|
# # ax[1, k].legend(loc='best')
|
||
|
# # ax[2, k].legend(loc='best')
|
||
|
# # ax[3, k].legend(loc='best')
|
||
|
# # # ax[0, k].plot(time, baseload[k], label='Base load')
|
||
|
# # # # ax[1, k].plot(time, Buyprice, label='Buyprice')
|
||
|
# # # ax[1, k].legend(loc='best')
|
||
|
# # # # ax[1, k].plot(time, Sellprice, label='Sellprice')
|
||
|
# # # # ax[1, k].legend(loc='best')
|
||
|
# # # # ax[2, k].plot(time, PEWH[k], label='EWH Power')
|
||
|
# # # # ax[2, k].plot(time, tettaEWHwat[k], label='EWH Temp')
|
||
|
# # # # ax[2, k].legend(loc='best')
|
||
|
# # # # ax[3, k].plot(time, PESSCharge[k],label='ESS Charging power')
|
||
|
# # # # ax[3, k].plot(time, PESSDischHome[k], label='ESS Disch to home')
|
||
|
# # # # ax[3, k].plot(time, PESSDischGrid[k], label='ESS Disch to grid')
|
||
|
# # # # ax[3, k].legend(loc='best')
|
||
|
# # # # ax[4, k].plot(time, PVHome[k], label='PV to Home')
|
||
|
# # # # ax[4, k].plot(time, PVGrid[k], label='PV to Grid')
|
||
|
# # # # ax[4, k].plot(time, PPV[k], label='All PV Production')
|
||
|
# # # # ax[4, k].legend(loc='best')
|
||
|
# # # # ax[5, k].plot(time, PEVCharge[k], label='EV Charging power')
|
||
|
# # # # ax[5, k].plot(time, PEVDischHome[k], label='EV Disch to home')
|
||
|
# # # # ax[5, k].plot(time, PEVDischGrid[k], label='EV Disch to grid')
|
||
|
# # # # ax[5, k].legend(loc='best')
|
||
|
# # # # ax[6, k].plot(time, EnergyESS[k], label='Energy of ESS')
|
||
|
# # # # ax[6, k].legend(loc='best')
|
||
|
# # # # ax[7, k].plot(time, EnergyEV[k], label='Energy of EV')
|
||
|
# # # # ax[7, k].legend(loc='best')
|
||
|
# # # # ax[4, k].set_xlabel('Time (step)')
|
||
|
# # # # ax[1, k].set_ylabel('Electricity price ($/W/h)')
|
||
|
# ax[0, 0].set_ylabel(
|
||
|
# 'Buy-Power1 (kW)') # showing only on the side of the graph if you want to put it on all need to use [0,k]
|
||
|
# ax[1, 0].set_ylabel('Buy-Power2 (kW)')
|
||
|
# ax[2, 0].set_ylabel('Sell-Power1 (kW)')
|
||
|
# ax[3, 0].set_ylabel('Sell-Power2 (kW)')
|
||
|
#
|
||
|
# ax[0, k].set_title(f"Home {k + 1}", fontsize='xx-small')
|
||
|
# fig.suptitle('Cost and Energy minimization', fontsize=16)
|
||
|
#
|
||
|
# plt.show()
|
||
|
# #
|
||
|
# # # if you want to plot each home alone in the figure
|
||
|
# # # for k in range(len(m2.h)):
|
||
|
# # # fig, ax = plt.subplots(2)
|
||
|
# # # ax[0].plot(time, pbuy[k], 'b-', label='Base')
|
||
|
# # # ax[0].plot(time, psell[k], 'g--', label='Controllable Load')
|
||
|
# # # ax[1].plot(time, Buyprice, label='buying price')
|
||
|
# # # ax[1].plot(time, Sellprice, 'y-+', label='selling price')
|
||
|
# # # ax[1].legend(loc='upper left')
|
||
|
# # # ax[1].set_xlabel('Time (hour)')
|
||
|
# # # ax[1].set_ylabel('Electricity price ($/W/h)')
|
||
|
# # # plt.suptitle(f"Home {k+1}")
|
||
|
# # # pl.tight_layout()
|
||
|
# # # plt.show()
|
||
|
# #
|
||
|
# # print('this is the end of code')
|
||
|
# # </editor-fold>
|
||
|
#
|
||
|
# # </editor-fold>
|
||
|
# </editor-fold>
|
||
|
|
||
|
# </editor-fold>
|
||
|
|
||
|
|
||
|
# Plot stacked bar chart
|
||
|
# fig11, ax = plt.subplots()
|
||
|
# labels = ('EV_Charging', 'EWH', 'Baseload', 'ESS_Charging')
|
||
|
# data_list = [PEVCharge2[1], PEWH2[1], baseload2[1], PESSChargeGrid2[1]]
|
||
|
# bottom = np.zeros(len(time))
|
||
|
# for label, data in zip(labels, data_list):
|
||
|
# ax.bar(time, data, bottom=bottom, label=label)
|
||
|
# bottom += data
|
||
|
# # Set plot labels and title
|
||
|
# ax.set_ylabel('Power (kW)')
|
||
|
# ax.set_xlabel('Time Step')
|
||
|
# ax.set_title('House Appliances (HEMS)_Energy Min')
|
||
|
# ax.legend(loc='upper left')
|
||
|
#
|
||
|
# # Show plot
|
||
|
# plt.show()
|
||
|
#
|
||
|
#
|
||
|
# fig111_en, ax = plt.subplots()
|
||
|
# labels_en = ('EV-Discharging_Grid','EV-Discharging_Home', 'PV_Grid', 'PV_Home','PV_Battery','ESS-Discharging_Grid','ESS-Discharging_Home')
|
||
|
# data_list_en = [PEVDischGrid[1],PEVDischHome[1], PVGrid[1],PVHome[1],PVBattery[1], PESSDischGrid[1],PESSDischHome[1]]
|
||
|
# bottom_en = np.zeros(len(time))
|
||
|
# for label_en, data_en in zip(labels_en, data_list_en):
|
||
|
# ax.bar(time, data_en, bottom=bottom_en, label=label_en)
|
||
|
# bottom2+= data_en
|
||
|
# # Set plot labels and title
|
||
|
# ax.set_ylabel('Power (kW)')
|
||
|
# ax.set_xlabel('Time Step')
|
||
|
# ax.set_title('Power of House Appliances export to Home/Grid (HEMS)-Energ Min')
|
||
|
# ax.legend(loc='upper left')
|
||
|
#
|
||
|
# # Show plot
|
||
|
# plt.show()
|
||
|
|
||
|
writer = pd.ExcelWriter('Result_HEMS.xlsx', engine='xlsxwriter')
|
||
|
time_excel = pd.DataFrame({'Time_Step': time})
|
||
|
time_excel.to_excel(writer, sheet_name='Power', startcol=0, index=False, header=True)
|
||
|
for h in range(len(m1.h)):
|
||
|
pload_excel = pd.DataFrame({f"Power_Demand {h + 1}": pl1_pbuy[h]})
|
||
|
pload_excel.to_excel(writer, sheet_name='Power', startcol=h + 1, index=False, header=True)
|
||
|
|
||
|
pload_excel = pd.DataFrame({f"Power_Sell {h + 1}": pl1_psell[h]})
|
||
|
pload_excel.to_excel(writer, sheet_name='Power', startcol=h + 4, index=False, header=True)
|
||
|
writer.save()
|