Upload files to "/"
This commit is contained in:
commit
30bd12534a
228
transformer_model.txt
Normal file
228
transformer_model.txt
Normal file
|
|
@ -0,0 +1,228 @@
|
|||
m3 = py.ConcreteModel()
|
||||
|
||||
# tranformer parameters
|
||||
m3.Transformer_Limit = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_data.Transformer_limit)),
|
||||
doc='Maximum power limit on transformer')
|
||||
m3.tran_life = py.Param(mutable=True, doc='Transformer life expectancy')
|
||||
m3.HST = py.Var(m1.t, bounds=(0, None), doc="windign hottest spot temp")
|
||||
m3.LOL = py.Var(m1.t, bounds=(0, None), doc="lost of life")
|
||||
# m3.FAA = py.Var(m1.t, bounds=(0, None), doc="HST")
|
||||
m3.DTOT = py.Var(m1.t, bounds=(0, None), doc="top oil rise temp")
|
||||
m3.DHST = py.Var(m1.t, bounds=(0, None), doc="winding HST rise over the top")
|
||||
m3.DTOTU = py.Var(m1.t, bounds=(0, None), doc="ultimate top-oil rise over the ambient temperature")
|
||||
m3.DHSTU = py.Var(m1.t, bounds=(0, None), doc="ultimate top-oil rise over the ambient")
|
||||
m3.Power_P_T = py.Var(m1.t, bounds=(0, None), doc="power from home to Transformer")
|
||||
m3.Power_T_P = py.Var(m1.t, bounds=(0, None), doc="power from transformer")
|
||||
m3.S_P_T = py.Var(m1.h, m1.t, within=py.Binary)
|
||||
m3.FAA = py.Var(m1.t, bounds=(0, None)) # The output of piecewise. Equivalent to model.y=piecewise(model.x)
|
||||
m3.n = py.Param(initialize=0.9,
|
||||
doc='top oil rise constant for transformer ') # 1, 0.9=ONOF tranformer, 0.8=ONON Transformer, it was 1
|
||||
m3.m = py.Param(initialize=1.6,
|
||||
doc='winding constant for transformer') # 1, 0.8 , 0.8 so for 2m will 2, 1.6 , 1.6 it was 2
|
||||
# <editor-fold desc="transformer parameter">
|
||||
m3.g1 = py.Param(initialize=0.1535,
|
||||
doc='time duration ') # i t was 0.4866, but 0.1535 if we consider delta t to 15 min if we consider 0.25 detla t then 0.0027
|
||||
# m3.g1 = py.Param(initialize=0.0027,doc='time duration ') # i t was 0.4866, if we consider 0.25 detla t then 0.0027
|
||||
m3.g2 = py.Param(initialize=0.8826,
|
||||
doc='time duration ') # it was 0.9998 but 0.8826 if we consider delta t to 15 min if we consider 0.25 detla t then 0.038
|
||||
# m3.g2 = py.Param(initialize=0.038, doc='time duration ') # it was 0.9998 but if we consider 0.25 detla t then 0.038
|
||||
m3.TR = py.Param(initialize=6, doc='Resistance ')
|
||||
m3.DTOR = py.Param(initialize=56, doc='e top-oil rise over ambient at the rated load ') # it was 56
|
||||
m3.DHR = py.Param(initialize=80, doc='hottest-spot rise over top-oil at the rated load ') # 24 and 80 it was 24
|
||||
m3.HSTL = py.Param(initialize=110, doc='time duration ')
|
||||
m3.t_TO = py.Param(initialize=180, doc=' in min ') # it was 90min in one paper and other paper 180 min
|
||||
m3.t_w = py.Param(initialize=48, doc=' in min ') # it was 7min in one paper and other paper 48 min
|
||||
m3.Ambient_Transformer = py.Param(m1.t, initialize=dict(zip(sheet_data.time, sheet_tr.celsius1)), mutable=True,
|
||||
doc='Buying Price')
|
||||
m3.P_Buy_extra = py.Param(m1.h, m1.t, initialize=0)
|
||||
m3.P_sell_extra = py.Param(m1.h, m1.t, initialize=0)
|
||||
delta_t = 15
|
||||
|
||||
|
||||
# <editor-fold desc="Transformer thermal Modeling">
|
||||
def HST_A(m3, i):
|
||||
return m3.HST[i] == m3.Ambient_Transformer[i] + m3.DTOT[i] + m3.DHST[i]
|
||||
|
||||
|
||||
m3.Const_Tr_a = py.Constraint(m1.t, rule=HST_A, doc='Hotest spot temperature')
|
||||
|
||||
|
||||
def HST_B(m3, i):
|
||||
# return m3.DTOT[i] == m3.DTOTU[i] * m3.g1
|
||||
return m3.DTOT[i] == m3.DTOTU[i] * (1 - py.exp(-delta_t / m3.t_TO.value))
|
||||
|
||||
|
||||
m3.Const_Tr_b = py.Constraint([m1.t_first.value], rule=HST_B, doc='Hotest spot temperature')
|
||||
|
||||
|
||||
def HST_C(m3, i):
|
||||
return m3.DTOT[i] == (m3.DTOTU[i] - m3.DTOT[i - 1]) * (1 - py.exp(-delta_t / m3.t_TO.value)) + m3.DTOT[i - 1]
|
||||
# return m3.DTOT[i] == (m3.DTOTU[i] - m3.DTOT[i - 1]) * m3.g1 + m3.DTOT[i - 1]
|
||||
|
||||
|
||||
m3.Const_Tr_c = py.Constraint(m1.t_second, rule=HST_C, doc='Hotest spot temperature')
|
||||
|
||||
|
||||
def HST_D(m3, i):
|
||||
return m3.DHST[i] == m3.DHSTU[i] * (1 - py.exp(-delta_t / m3.t_w.value))
|
||||
# return m3.DHST[i] == m3.DHSTU[i] * m3.g2.value
|
||||
|
||||
|
||||
m3.Const_Tr_d = py.Constraint([m1.t_first.value], rule=HST_D, doc='Hotest spot temperature')
|
||||
|
||||
|
||||
def HST_E(m3, i):
|
||||
# return m3.DHST[i] == (m3.DHSTU[i] - m3.DHST[i - 1]) * m3.g2.value + m3.DHST[i - 1]
|
||||
return m3.DHST[i] == (m3.DHSTU[i] - m3.DHST[i - 1]) * (1 - py.exp(-delta_t / m3.t_w.value)) + m3.DHST[i - 1]
|
||||
|
||||
|
||||
m3.Const_Tr_e = py.Constraint(m1.t_second, rule=HST_E, doc='Hotest spot temperature')
|
||||
|
||||
|
||||
# m = 0.8, n =0.9 or .8 ,.8, or 1,1
|
||||
def HST_F(m3, i):
|
||||
# return m3.DTOTU[i] == (
|
||||
# m3.DTOR.value * ((m3.Power_P_T[i].value ** 2) * m3.TR.value + 1)) / (
|
||||
# m3.TR.value + 1)
|
||||
return m3.DTOTU[i] == m3.DTOR.value * (
|
||||
(((m1.Power_Total_power[i].value / m3.Transformer_Limit[i]) ** 2) * m3.TR.value + 1) / (
|
||||
m3.TR.value + 1)) ** m3.n.value
|
||||
|
||||
|
||||
m3.Const_Tr_f = py.Constraint(m1.t, rule=HST_F, doc='Hotest spot temperature')
|
||||
|
||||
|
||||
def HST_G(m3, i):
|
||||
return m3.DHSTU[i] == m3.DHR.value * (
|
||||
((m1.Power_Total_power[i].value) / m3.Transformer_Limit[i]) ** m3.m.value) # was 2 , 2*0.8 = 1.6
|
||||
|
||||
|
||||
m3.Const_Tr_g = py.Constraint(m1.t, rule=HST_G, doc='Hotest spot temperature')
|
||||
|
||||
|
||||
# def HSTL(m3, i):
|
||||
# return m3.HST[i] <= m3.HSTL.value
|
||||
#
|
||||
#
|
||||
# m3.Const_Tr_i = py.Constraint(m1.t, rule=HSTL, doc='Hotest upeer limit')
|
||||
#
|
||||
|
||||
def Lol1(m3, i):
|
||||
return m3.FAA[i] == py.exp(15000 / 383 - (15000 / (m3.HST[i] + 273)))
|
||||
|
||||
|
||||
m3.Const_Tr_j = py.Constraint(m1.t, rule=Lol1, doc='Aging factor')
|
||||
|
||||
|
||||
def Lol(m3, i):
|
||||
return m3.LOL[i] == (sum(m3.FAA[y] for y in m1.t) * m1.time_d.value) / 180000 # it was m1.time_d.value
|
||||
# return m3.LOL[i] == (m3.FAA[i] * m1.time_d.value) / 180000
|
||||
|
||||
|
||||
m3.Const_Tr_k = py.Constraint(m1.t, rule=Lol, doc='LOST OF LIFE transformer')
|
||||
|
||||
|
||||
# LDRv(t)=E= SUM(td, LDR(t,td)) ;
|
||||
# MP(t).. Ptot(t)=E= P(t)-LLC(t)-b*PNSS-SUM(td, LDR(t,td))+SUM(td, LDR(td,t)) ;
|
||||
# C_HST(t).. HST(t)=l=HSTL;# GAMS code
|
||||
# C_LC(t).. LLC(t)=l=Ptot(t)-LLC(t);# HST1(t).. HST(t)=e=AT(t)+DTOT(t)++(t);
|
||||
|
||||
|
||||
# OF .. Cost =E= SUM(t, LLC(t))*CLC + b*CSW + SUM(t, LDRv(t))*CDR ;
|
||||
|
||||
|
||||
def objective_rule(m3):
|
||||
return m3.LOL[i]
|
||||
# return m3.LOL[i]
|
||||
|
||||
|
||||
m3.obj_HEMS = py.Objective(rule=objective_rule, sense=py.minimize, doc='Definition of objective function')
|
||||
# <editor-fold desc="Solver">
|
||||
# m3.write('mtr12.lp', io_options={'symbolic_solver_labels': True})
|
||||
opt = py.SolverFactory('ipopt')
|
||||
# opt = py.SolverFactory('mindtpy').solve(m3, mip_solver='cplex', nlp_solver='ipopt', tee=True)
|
||||
# opt = py.SolverFactory('couenne ')
|
||||
# opt.options["mipgap"] = 1 # 0.155 for load and 0.8 for other load new one: 0.045
|
||||
result = opt.solve(m3, tee=True)
|
||||
# result = opt.solve(m2)
|
||||
|
||||
# this is for calling the neos
|
||||
# os.environ['NEOS_EMAIL'] = 'sadamengr15@gmail.com'
|
||||
# opt = py.SolverManagerFactory('neos')
|
||||
# result = opt.solve(m2, opt='cplex')
|
||||
|
||||
# m2.pprint()
|
||||
# print(result)
|
||||
|
||||
# print( "results.solver.status =", result.solver.status)
|
||||
# print("results.solver.termination_condition =", result.solver.termination_condition )
|
||||
# print("SolverStatus.warning =", result.warning)
|
||||
print(result)
|
||||
|
||||
|
||||
|
||||
# # # <editor-fold desc="Ploting">
|
||||
# #
|
||||
# # # Plotting for automatic ploting
|
||||
time = [i for i in m1.t]
|
||||
|
||||
HST = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
DHST = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
DTOT = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
Transformer_Limit = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
Power_P_T = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
Power_T_P = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
Load_Transformer = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
Ambient_Transformer = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
FAA = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
LOL = [] # cost = [i for i in m2.c.values()] <<< check this
|
||||
FAA_1 = [1] * 96
|
||||
HST_ref = [110] * 96
|
||||
for i in m3.LOL:
|
||||
LOL.append(py.value(m3.LOL[i]))
|
||||
for i in m3.FAA:
|
||||
FAA.append(py.value(m3.FAA[i]))
|
||||
for i in m3.Ambient_Transformer:
|
||||
Ambient_Transformer.append(py.value(m3.Ambient_Transformer[i]))
|
||||
for i in m3.HST:
|
||||
HST.append(py.value(m3.HST[i]))
|
||||
for i in m3.DTOT:
|
||||
DTOT.append(py.value(m3.DTOT[i]))
|
||||
for i in m3.DHST:
|
||||
DHST.append(py.value(m3.DHST[i]))
|
||||
# for i in m1.Power_P_T:
|
||||
# Power_P_T.append(py.value(m1.Power_P_T[i]))
|
||||
for i in m3.Transformer_Limit:
|
||||
Transformer_Limit.append(py.value(m3.Transformer_Limit[i]))
|
||||
for i in m1.Power_Total_power:
|
||||
Load_Transformer.append(py.value(m1.Power_Total_power[i]))
|
||||
fig_T1, ax = plt.subplots(3, constrained_layout=True, ) # , sharey='row'
|
||||
gs = fig_T1.add_gridspec(hspace=0, wspace=0)
|
||||
# ax[0].plot(time, Transformer_Limit, color="red")
|
||||
# ax[0].bar(time, Power_P_T, color="blue")
|
||||
# ax[0].legend(['Tranformer_limit'])
|
||||
ax[0].plot(time, Transformer_Limit, color="red",label="Transformer limit")
|
||||
ax[0].bar(time, Load_Transformer, color="blue",label="Load on Transformer")
|
||||
ax[0].legend()
|
||||
ax[1].plot(time, HST, label='HST', color="red")
|
||||
# ax[1].plot(time, DHST, label='DHST', color="blue")
|
||||
# ax[1].plot(time, DTOT, label='DTOT', color="green")
|
||||
ax[1].plot(time, HST_ref, label='HST_ref',linestyle=':', color="black")
|
||||
ax[1].legend()
|
||||
ax[0].set_ylabel('Load_Transformer(kW)')
|
||||
ax[2].plot(time, FAA, label='FAA', color="red")
|
||||
ax[2].plot(time, FAA_1, label='FAA=1',linestyle=':' ,color="black")
|
||||
ax[2].legend()
|
||||
# ax[3].plot(time, LOL, label='% LOL', color="blue")
|
||||
# ax[4].bar(HST, FAA, label='LOL', color="blue")
|
||||
# ax[5].bar(HST, LOL, label='LOL', color="blue")
|
||||
ax[0].set_ylabel('Load_Transformer (kW)')
|
||||
ax[1].set_ylabel('Temp')
|
||||
ax[2].set_ylabel('Acc aging factor')
|
||||
# ax[3].set_ylabel('LoL %')
|
||||
# ax[4].set_ylabel('HST Vs FAA')
|
||||
# ax[5].set_ylabel('HST Vs %LoL')
|
||||
fig_T1.suptitle('Transformer connected to HEMS', fontsize=16)
|
||||
# #
|
||||
# pl.tight_layout()
|
||||
plt.show()
|
||||
Loading…
Reference in New Issue
Block a user