In [3]:
from __future__ import division
from pulp import *
import os
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.dates as mdates
import matplotlib.mlab as mlab
import pandas as pd
import csv
import datetime
import sqlite3
import math
import glob
import seaborn as sns
%matplotlib inline
from pandas import *
from scipy import optimize
from scipy import stats
from datetime import date, timedelta, datetime
import timeit
from collections import OrderedDict
from IPython.display import clear_output
from IPython.display import HTML
from matplotlib.ticker import MultipleLocator
from shutil import copyfile
In [13]:
# df_header = pd.read_csv('results/FlexDrivers/NewEVStatus/2015-07-01.csv')
filenames = glob.glob('results/FlexDrivers/NewEVStatus/alcopark/EventSessions/'+'*.csv')
dfs = pd.DataFrame()
for filename in filenames:
df = pd.read_csv(filename)
df.index = df[df.columns[0]]
dfs = pd.concat([dfs,df])
dfs.drop(dfs.columns[[0]], axis=1, inplace=True)
#dfs.to_csv('results/FlexDrivers/NewEVStatus/alcopark/' + '2015-EV-Sessions.csv')
In [14]:
station_sessions = []
station_sessions = dfs.copy()
In [15]:
def query_EVSE_eventids(startdate, enddate):
# Query each event, for each event run the optimization based on corresponding energy cost to minimize
# the sum of energy cost every 15 minutes.
print('Start to query EV power at "Alcopark"')
# Get distinct event_ids for
conn = sqlite3.connect("sessionData.db")
cursor = conn.cursor()
query=("SELECT DISTINCT(EVENT_ID) FROM sessiondata")
query_cond1=("where ENERGY>0 and ENERGY!='' AND")
query_cond2=("STATION LIKE '%ALCO%' AND") # ALCOPARK
# PARKING AT NIGHT, FLEET VEHICLES
# query_cond3=("(((CAST(strftime('%j',SESSION_END_TIME) as INT) - CAST(strftime('%j',SESSION_START_TIME) as INT)) > 0) AND")
query_cond4=("(strftime('%Y-%m-%d',SESSION_START_TIME) between ") # MONTH IN JAN
cur = cursor.execute(query+" "+query_cond1+" "+query_cond2+" "+query_cond4+"'"+startdate+"'"+' AND '+"'"+enddate+"'"+')')
rows = cur.fetchall()
event_ids=[]
for row in rows: event_ids.append(row[0])
event_ids=np.array(event_ids)
conn.close()
# query EVSE data from SQL
EVSE_Power = pd.DataFrame()
index = pd.date_range(startdate, enddate, freq='15min')
EVSE_Power = pd.DataFrame(index=index)
tcts=len(event_ids)
for event_id in event_ids:
# print('This is the current EVENT_ID {}'.format(event_id))
# if tcts % 30 ==0: clear_output()
# start = timeit.default_timer()
conn = sqlite3.connect("minData.db")
cursor = conn.cursor()
# FORM A QUERY TO GET RELEVANT DATA, MAKE SURE TO FLOOR THE STARTING OR ENDING POINT
# THIS WUERY WILL RETURN INTERVAL_ENERGY IN KWH, THIS REALLY IS WHAT WE ARE INTERESTED IN FOR OPTIMIZATION
# WE WILL NEED THE TIME FIELD TO GET THE CORRESPONDING VALUES OF ENERGY COST
# THE REST WILL BE USED WHEN POPULATING THE MINDATA DATABASE
query=("SELECT INTERVAL_START_DATE AS START,INTERVAL_END_DATE AS END,AVERAGE_POWER,ENERGY FROM minData WHERE ENERGY!='' AND EVENT_ID=")
query_cond=str(event_id)
cur = cursor.execute(query+query_cond)
df = pd.DataFrame(cur.fetchall())
df.columns=['start','end','power','cumulative_energy']
df['StartDate'] = df['start'].apply(lambda x: datetime.strptime(x,'%Y-%m-%d %H:%M:%S'))
df.index = df['StartDate']
try:
df = df.resample('15min').bfill()
df = df.rename(columns = {'power':event_id})
EVSE_Power = pd.concat([EVSE_Power,pd.DataFrame(df[event_id])], axis = 1)
except ValueError:
print('Oops! That was no valid Event ID: ' + str(event_id))
except KeyError:
print ('Void Dict no valid Event ID: ' + str(event_id))
return EVSE_Power
In [16]:
# consider the holidays
holidays = pd.read_csv('holiday.csv', usecols=[0])
holidays['date'] = pd.to_datetime(pd.Series(holidays['date']), format='%m/%d/%y')
holidays['day'] = holidays.date.apply(lambda x: x.strftime('%Y-%m-%d'))
# utility tariff, E19-TOU
utility = []
utility = pd.read_csv('PGE_E19.csv', usecols=[1])
# meter charge
meter_charge = float("{0:.2f}".format(utility['Secondary'][0]))
# demand charge
summer_onpeak_demand = float("{0:.5f}".format(utility['Secondary'][1]))
summer_midpeak_demand = float("{0:.5f}".format(utility['Secondary'][2]))
summer_monthly_demand = float("{0:.5f}".format(utility['Secondary'][3]))
winter_midpeak_demand = float("{0:.5f}".format(utility['Secondary'][4]))
winter_monthly_demand = float("{0:.5f}".format(utility['Secondary'][5]))
# energy charge
summer_onpeak_energy = float("{0:.5f}".format(utility['Secondary'][6]))
summer_midpeak_energy = float("{0:.5f}".format(utility['Secondary'][7]))
summer_offpeak_energy = float("{0:.5f}".format(utility['Secondary'][8]))
winter_midpeak_energy = float("{0:.5f}".format(utility['Secondary'][9]))
winter_offpeak_energy = float("{0:.5f}".format(utility['Secondary'][10]))
# PDP charge and credits
cpp_charge = utility['Secondary'][11]
credits_onpeak_demand = utility['Secondary'][12]
credits_midpeak_demand = utility['Secondary'][13]
# credits_onpeak_energy = utility['Secondary'][14]
# credits_midpeak_energy = utility['Secondary'][15]
# Get cost value of energy charges for time
def getCostValueForTime(month,time,weekday):
if (month in [1,2,3,4,11,12]) & (weekday > 0) & (weekday < 6):
if (time >= '08:30') & (time < '21:30'):
return winter_midpeak_energy
elif (time >= '00:00') & (time < '08:30'):
return winter_offpeak_energy
elif time >= '21:30':
return 0
elif (month in [1,2,3,4,11,12]) & (weekday in [6,7]):
return winter_offpeak_energy
elif (month in [5,6,7,8,9,10]) & (weekday > 0) & (weekday < 6):
if (time >= '00:00') & (time < '08:30'):
return summer_offpeak_energy
elif (time >= '08:30') & (time < '12:00'):
return summer_midpeak_energy
elif (time >= '12:00') & (time < '18:00'):
return summer_onpeak_energy
elif (time >= '18:00') & (time < '21:30'):
return summer_midpeak_energy
elif (time >= '21:30') & (time <= '23:45'):
return summer_offpeak_energy
elif (month in [5,6,7,8,9,10]) & (weekday in [6,7]):
return summer_offpeak_energy
else:
print "Error in month, time and weekday"
In [17]:
ts2015 = pd.read_csv('results/FlexDrivers/NewEVStatus/2015TimeStamp.csv')
ts2015.index = ts2015.ix[:,0].apply(lambda x: datetime.strptime(x,'%m/%d/%y %H:%M'))
ts2015.drop(ts2015.columns[[0]], axis=1, inplace=True)
ts2015['day'] = map(lambda x: x.strftime('%Y-%m-%d'),ts2015.index)
ts2015['month'] = map(lambda x: x.strftime('%m'),ts2015.index)
ts2015['hour'] = map(lambda x: x.strftime('%H'),ts2015.index)
ts2015['hour'] = ts2015.hour.astype(int)
ts2015['month'] = ts2015.month.astype(int)
In [18]:
import random
def RandomEVOpt(station_sessions, driver_pct, event_pct):
NonFleetEV_Sessions = station_sessions.loc[station_sessions.driver_id.apply(lambda x: (x not in [146547,294291])),:]
grp_df2 = pd.DataFrame(NonFleetEV_Sessions.groupby(['driver_id'])['event_id'].size())
grp_df2.columns = ['session_num']
FleetEV_ID = list(station_sessions.loc[station_sessions.driver_id.apply(lambda x: (x in [146547,294291])),'event_id'].values)
NonFleetEV_ID = list(station_sessions.loc[station_sessions.driver_id.apply(lambda x: (x not in [146547,294291])),'event_id'].values)
NonFleetEV_50 = list(grp_df2.sort_values(['session_num'],ascending=False).head(n=int(len(grp_df2)*driver_pct)).index)
NonFleetEV_ID_50 = list(station_sessions.loc[station_sessions.driver_id.apply(lambda x: (x in NonFleetEV_50)),'event_id'].values)
NonFleetEV_ID_50_Random = random.sample(NonFleetEV_ID_50,int(len(NonFleetEV_ID_50)*event_pct))
EV_ID_Random = FleetEV_ID + NonFleetEV_ID_50_Random
return EV_ID_Random
EV_ID_Random = RandomEVOpt(station_sessions,1.0,1.0)
In [19]:
len(EV_ID_Random)
Out[19]:
In [20]:
def hasSlack(eventid,DailyPower):
DayPower = DailyPower.fillna(0)
oldpower = DailyPower.loc[DailyPower[eventid] >= 0,[eventid,'ts']]
newpower = DayPower.loc[DayPower[eventid] >= 0,[eventid]]
power2 = newpower[eventid].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS_ref = oldpower['ts'][0]
RS_ref = oldpower['ts'][-1]
if M == 0:
return False
elif M >= (RS_ref - LS_ref + 1):
return False
else:
return True
def flexEvent(eventid, DailyPower, past_event):
flex_event = []
DayPower = DailyPower.fillna(0)
oldpower = DailyPower.loc[DailyPower[eventid] >= 0,[eventid,'ts']]
newpower = DayPower.loc[DayPower[eventid] >= 0,[eventid]]
power2 = newpower[eventid].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS_ref = oldpower['ts'][0]
RS_ref = oldpower['ts'][-1]
past_event.append(eventid)
for event in past_event:
oldpower = DailyPower.loc[(DailyPower['ts'] >= LS_ref) & (DailyPower[event] >= 0),[event,'ts']]
newpower = DayPower.loc[(DailyPower['ts'] >= LS_ref) & (DailyPower[event] >= 0),[event]]
power2 = newpower[event].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
if len(oldpower) > 0:
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
if M == 0:
pass
elif M >= (RS - LS + 1):
pass
else:
flex_event.append(event)
else:
pass
return flex_event, LS_ref, RS_ref
In [10]:
def OptimizeCostV5(DailyPower, LS_ref, eventid, past_event, opt_event, onpeak, midpeak, peak):
DailyPower2 = DailyPower.copy()
DailyPower2.loc[DailyPower2.index[0:LS_ref],opt_event] = np.nan
DayPower = DailyPower2.fillna(0)
# opt_event.append('power')
# Create matrix of slots
X = LpVariable.dicts("ChargingStatus", [(k, i, j) for k in range(len(opt_event)) for i in range(96) for j in range(96)], 0, 1, LpBinary)
costs=[]
for tim in DayPower.index:
costs.append(getCostValueForTime(tim.month,tim.strftime('%H:%M'),tim.isoweekday()))
cost={}
for k in xrange(len(opt_event)):
for i in xrange(96):
for j in xrange(96):
cost[k,i,j]=costs[j]
# Create model objective
model = pulp.LpProblem("SmartCharging",pulp.LpMinimize)
# Create model variables
min_curMidPeak = LpVariable("Min MidPeak",0)
min_curOnPeak = LpVariable("Min OnPeak",0)
min_curPeak = LpVariable("Min Peak",0)
tim = DayPower.index[0]
# Create model constrains
for k in xrange(len(opt_event)):
oldpower = DailyPower2.loc[DailyPower2[opt_event[k]] >= 0,[opt_event[k],'ts']]
newpower = DayPower.loc[DayPower[opt_event[k]] >= 0,[opt_event[k]]]
power2 = newpower[opt_event[k]].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
# Set the precedence and the session duration constraints
for i in xrange(LS):
model += (lpSum(X[k,i,j] for j in range(96)) == 0)
# Down zeros
for i in xrange(LS+M,96):
model += (lpSum(X[k,i,j] for j in range(96)) == 0)
# Left zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(LS)) == 0)
# Right zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(RS+1,96)) == 0)
# Add the constraints
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(LS,RS+1)) == 1)
for i in xrange(LS,LS+M-1):
# model += (lpSum(X[i,j]*(j+1) for j in range(LS,RS+1)) == lpSum(X[i+1,j]*(j+1) for j in range(LS,RS+1))-1)
model += (lpSum(X[k,i,j]*(j+1) for j in range(LS,RS+1)) <= lpSum(X[k,i+1,j]*(j+1) for j in range(LS,RS+1))-1) # discrete 15-min charging load blocks
# Re-ordering charging slot for next round optimization
DayPower.loc[DayPower.index[0:96],opt_event[k]] = [0] * LS + positive_power + [0] * (96-LS-M)
# Create model objective
if (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [1,2,3,4,5]):
model += lpSum((DayPower[opt_event[k]][i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak * 17.33 + min_curMidPeak*0.13
for j in xrange(34,86):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curMidPeak
for j in xrange(96):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curPeak
elif (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [6,7]):
model += lpSum((DayPower[opt_event[k]][i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak * 17.33
for j in xrange(96):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curPeak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [1,2,3,4,5]):
model += lpSum((DayPower[opt_event[k]][i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curOnPeak * 18.74 + min_curPeak * 17.33 + min_curMidPeak * 5.23
for j in (range(34,48)+range(72,86)):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curMidPeak
for j in xrange(96):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curPeak
for j in xrange(48,72):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curOnPeak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [6,7]):
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33
for j in xrange(96):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curPeak
else:
print "Error in month, day and time"
# model.solve()
model.solve(PULP_CBC_CMD(fracGap = 0.05, maxSeconds = 300, threads = None))
if model.status == 1:
for k in xrange(len(opt_event)):
oldpower = DailyPower2.loc[DailyPower2[opt_event[k]] >= 0,[opt_event[k],'ts']]
newpower = DayPower.loc[DayPower[opt_event[k]] >= 0,[opt_event[k]]]
power2 = newpower[opt_event[k]].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
slots = xrange(1,len(power2)+1)
sol=[]
for j in xrange(96):
sol.append(sum(X[k,i,j].varValue*(j) for i in xrange(LS,LS+M)))
# corSlots=[slots[int(s)] for s in sol]
ct=0
power=[p*0 for p in power2]
for s in sol:
if s>0:
power[int(s)]=positive_power[ct]
ct=ct+1
else:
continue
DailyPower.loc[DailyPower.index[LS:RS+1],opt_event[k]] = power[LS:RS+1]
else:
print('The solution is suboptimal')
return DailyPower
In [11]:
# run the optimization for each month's EV charging sessions
for month in range(7,8):
# monthly date range
StartDate = ts2015[ts2015['month'] == month].day[0]
EndDate = ts2015[ts2015['month'] == month].day[-1]
rng = map(lambda x: x.strftime('%Y-%m-%d'),pd.date_range(StartDate, EndDate, freq='D'))
EndDate2 = (datetime.strptime(EndDate,'%Y-%m-%d')+timedelta(days=1)).strftime('%Y-%m-%d')
# EVSE_Power = query_EVSE_eventids(StartDate,EndDate2)
# EVSE_Power = EVSE_Power[StartDate:EndDate]
EVSE_Power = pd.read_csv('results/FlexDrivers/NewEVStatus/2015EV/'+'2015-'+str(month)+'Original.csv')
EVSE_Power.index = EVSE_Power.ix[:,0].apply(lambda x: datetime.strptime(x,'%Y-%m-%d %H:%M:%S'))
EVSE_Power.drop(EVSE_Power.columns[[0]], axis=1, inplace=True)
# initialize the past peak demand values
DailyPower = []
peak=0
onpeak=0
midpeak=0
monthpeak = []
monthonpeak = []
monthmidpeak = []
directory = 'results/FlexDrivers/NewEVStatus/OptV5/' + str(month)
if not os.path.exists(directory):
os.makedirs(directory)
# run the daily optimization
for day in rng[25:26]:
DailyPower = EVSE_Power[day]
DailyPower = DailyPower.replace('',np.nan)
DailyPower = DailyPower.dropna(axis=1,how='all')
DriverOnList = DailyPower.columns
N = len(DriverOnList)
opt_event = []
past_event = []
NonSlack_event = []
DailyPower['power'] = 0
DailyPower['NonControlledPower'] = DailyPower.fillna(0).sum(axis=1)
DailyPower['ts'] = range(96) # Add the tag of time step for each event id
if DailyPower.empty:
print 'EV Data on ' + day + ' is empty!'
DailyPower.to_csv(directory + '/' + day +'.csv')
else:
for eventid in DriverOnList[::-1]:
if (int(eventid) in EV_ID_Random) & (hasSlack(eventid,DailyPower)):
opt_event, LS_ref, RS_ref = flexEvent(eventid, DailyPower, past_event)
df2 = DailyPower[opt_event].fillna(0).sum(axis=1)
df2[LS_ref:96] = 0
DailyPower['power'] = 0
DailyPower['power'] = DailyPower[NonSlack_event].fillna(0).sum(axis=1) + \
DailyPower[[x for x in past_event if x not in opt_event]].fillna(0).sum(axis=1) + df2
DailyPower = OptimizeCostV5(DailyPower, LS_ref, eventid, past_event, opt_event, onpeak, midpeak, peak)
print 'Optimization in ' + str(eventid)
else:
NonSlack_event.append(eventid)
print str(eventid) + ' has no slack for charging load move!'
DailyPower['power'] = 0
DailyPower['power'] = DailyPower[NonSlack_event].fillna(0).sum(axis=1) + DailyPower[past_event].fillna(0).sum(axis=1)
print 'Optimization on ' + day + ' is done successfully'
DailyPower.to_csv(directory +'/'+ day +'-V5.csv')
# EVSE_Power[StartDate:EndDate].to_csv('results/FlexDrivers/NewEVStatus/OptV5/'+'2015-'+str(month)+'Original.csv')
In [21]:
def OptimizeCostV6(DailyPower, LS_ref, eventid, past_event, opt_event, onpeak, midpeak, peak):
DailyPower2 = DailyPower.copy()
DailyPower2.loc[DailyPower2.index[0:LS_ref],opt_event] = np.nan
DayPower = DailyPower2.fillna(0)
# opt_event.append('power')
# Create matrix of slots
X = LpVariable.dicts("ChargingStatus", [(k, i, j) for k in range(len(opt_event)) for i in range(96) for j in range(96)], 0, 1, LpBinary)
costs=[]
for tim in DayPower.index:
costs.append(getCostValueForTime(tim.month,tim.strftime('%H:%M'),tim.isoweekday()))
cost={}
for k in xrange(len(opt_event)):
for i in xrange(96):
for j in xrange(96):
cost[k,i,j]=costs[j]
# Create model objective
model = pulp.LpProblem("SmartCharging",pulp.LpMinimize)
# Create model variables
min_curMidPeak = LpVariable("Min MidPeak",0)
min_curOnPeak = LpVariable("Min OnPeak",0)
min_curPeak = LpVariable("Min Peak",0)
tim = DayPower.index[0]
# Create model constrains
for k in xrange(len(opt_event)):
oldpower = DailyPower2.loc[DailyPower2[opt_event[k]] >= 0,[opt_event[k],'ts']]
newpower = DayPower.loc[DayPower[opt_event[k]] >= 0,[opt_event[k]]]
power2 = newpower[opt_event[k]].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
# Set the precedence and the session duration constraints
for i in xrange(LS):
model += (lpSum(X[k,i,j] for j in range(96)) == 0)
# Down zeros
for i in xrange(LS+M,96):
model += (lpSum(X[k,i,j] for j in range(96)) == 0)
# Left zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(LS)) == 0)
# Right zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(RS+1,96)) == 0)
# Add the constraints
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(LS,RS+1)) == 1)
for i in xrange(LS,LS+M-1):
# model += (lpSum(X[i,j]*(j+1) for j in range(LS,RS+1)) == lpSum(X[i+1,j]*(j+1) for j in range(LS,RS+1))-1)
model += (lpSum(X[k,i,j]*(j+1) for j in range(LS,RS+1)) <= lpSum(X[k,i+1,j]*(j+1) for j in range(LS,RS+1))-1) # discrete 15-min charging load blocks
# Re-ordering charging slot for next round optimization
DayPower.loc[DayPower.index[0:96],opt_event[k]] = [0] * LS + positive_power + [0] * (96-LS-M)
# Create model objective
if (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [1,2,3,4,5]):
model += lpSum((DayPower[opt_event[k]][i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak * 17.33
for j in xrange(34,86):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curMidPeak
for j in xrange(96):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curPeak
elif (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [6,7]):
model += lpSum((DayPower[opt_event[k]][i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak * 17.33
for j in xrange(96):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curPeak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [1,2,3,4,5]):
model += lpSum((DayPower[opt_event[k]][i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curOnPeak * 18.74 + min_curPeak * 17.33 + min_curMidPeak * 5.23
for j in (range(34,48)+range(72,86)):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curMidPeak
for j in xrange(96):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curPeak
for j in xrange(48,72):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curOnPeak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [6,7]):
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33
for j in xrange(96):
model += lpSum(DayPower[opt_event[k]][i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'][j] <= min_curPeak
else:
print "Error in month, day and time"
# model.solve()
model.solve(PULP_CBC_CMD(fracGap = 0.05, maxSeconds = 300, threads = None))
if model.status == 1:
for k in xrange(len(opt_event)):
oldpower = DailyPower2.loc[DailyPower2[opt_event[k]] >= 0,[opt_event[k],'ts']]
newpower = DayPower.loc[DayPower[opt_event[k]] >= 0,[opt_event[k]]]
power2 = newpower[opt_event[k]].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
slots = xrange(1,len(power2)+1)
sol=[]
for j in xrange(96):
sol.append(sum(X[k,i,j].varValue*(j) for i in xrange(LS,LS+M)))
# corSlots=[slots[int(s)] for s in sol]
ct=0
power=[p*0 for p in power2]
for s in sol:
if s>0:
power[int(s)]=positive_power[ct]
ct=ct+1
else:
continue
DailyPower.loc[DailyPower.index[LS:RS+1],opt_event[k]] = power[LS:RS+1]
else:
print('The solution is suboptimal')
return DailyPower
In [ ]:
# run the optimization for each month's EV charging sessions
for month in range(5,11):
# monthly date range
StartDate = ts2015[ts2015['month'] == month].day[0]
EndDate = ts2015[ts2015['month'] == month].day[-1]
rng = map(lambda x: x.strftime('%Y-%m-%d'),pd.date_range(StartDate, EndDate, freq='D'))
EndDate2 = (datetime.strptime(EndDate,'%Y-%m-%d')+timedelta(days=1)).strftime('%Y-%m-%d')
# EVSE_Power = query_EVSE_eventids(StartDate,EndDate2)
# EVSE_Power = EVSE_Power[StartDate:EndDate]
EVSE_Power = pd.read_csv('results/FlexDrivers/NewEVStatus/2015EV/'+'2015-'+str(month)+'Original.csv')
EVSE_Power.index = EVSE_Power.ix[:,0].apply(lambda x: datetime.strptime(x,'%Y-%m-%d %H:%M:%S'))
EVSE_Power.drop(EVSE_Power.columns[[0]], axis=1, inplace=True)
# initialize the past peak demand values
DailyPower = []
peak=0
onpeak=0
midpeak=0
monthpeak = []
monthonpeak = []
monthmidpeak = []
directory = 'results/FlexDrivers/NewEVStatus/OptV6/' + str(month)
if not os.path.exists(directory):
os.makedirs(directory)
# run the daily optimization
for day in rng:
DailyPower = EVSE_Power[day]
DailyPower = DailyPower.replace('',np.nan)
DailyPower = DailyPower.dropna(axis=1,how='all')
DriverOnList = DailyPower.columns
N = len(DriverOnList)
opt_event = []
past_event = []
NonSlack_event = []
DailyPower['power'] = 0
DailyPower['NonControlledPower'] = DailyPower.fillna(0).sum(axis=1)
DailyPower['ts'] = range(96) # Add the tag of time step for each event id
if DailyPower.empty:
print 'EV Data on ' + day + ' is empty!'
DailyPower.to_csv(directory + '/' + day +'.csv')
else:
for eventid in DriverOnList[::-1]:
if (int(eventid) in EV_ID_Random) & (hasSlack(eventid,DailyPower)):
opt_event, LS_ref, RS_ref = flexEvent(eventid, DailyPower, past_event)
df2 = DailyPower[opt_event].fillna(0).sum(axis=1)
df2[LS_ref:96] = 0
DailyPower['power'] = 0
DailyPower['power'] = DailyPower[[x for x in past_event if x not in opt_event]].fillna(0).sum(axis=1) + df2
DailyPower = OptimizeCostV6(DailyPower, LS_ref, eventid, past_event, opt_event, onpeak, midpeak, peak)
print 'Optimization in ' + str(eventid)
else:
NonSlack_event.append(eventid)
opt_event, LS_ref, RS_ref = flexEvent(eventid, DailyPower, past_event)
df2 = DailyPower[opt_event].fillna(0).sum(axis=1)
df2[LS_ref:96] = 0
DailyPower['power'] = 0
DailyPower['power'] = DailyPower[[x for x in past_event if x not in opt_event]].fillna(0).sum(axis=1) + df2
DailyPower = OptimizeCostV6(DailyPower, LS_ref, eventid, past_event, opt_event, onpeak, midpeak, peak)
print str(eventid) + ' has no slack for charging load move!'
print 'Optimization in ' + str(eventid)
DailyPower['power'] = 0
DailyPower['power'] = DailyPower[past_event].fillna(0).sum(axis=1)
print 'Optimization on ' + day + ' is done successfully'
DailyPower.to_csv(directory +'/'+ day +'-V6.csv')
# EVSE_Power[StartDate:EndDate].to_csv('results/FlexDrivers/NewEVStatus/OptV5/'+'2015-'+str(month)+'Original.csv')
In [23]:
# run the optimization for each month's EV charging sessions
for month in (range(4,5)+range(11,13)):
# monthly date range
StartDate = ts2015[ts2015['month'] == month].day[0]
EndDate = ts2015[ts2015['month'] == month].day[-1]
rng = map(lambda x: x.strftime('%Y-%m-%d'),pd.date_range(StartDate, EndDate, freq='D'))
EndDate2 = (datetime.strptime(EndDate,'%Y-%m-%d')+timedelta(days=1)).strftime('%Y-%m-%d')
# EVSE_Power = query_EVSE_eventids(StartDate,EndDate2)
# EVSE_Power = EVSE_Power[StartDate:EndDate]
EVSE_Power = pd.read_csv('results/FlexDrivers/NewEVStatus/2015EV/'+'2015-'+str(month)+'Original.csv')
EVSE_Power.index = EVSE_Power.ix[:,0].apply(lambda x: datetime.strptime(x,'%Y-%m-%d %H:%M:%S'))
EVSE_Power.drop(EVSE_Power.columns[[0]], axis=1, inplace=True)
# initialize the past peak demand values
DailyPower = []
peak=0
onpeak=0
midpeak=0
monthpeak = []
monthonpeak = []
monthmidpeak = []
directory = 'results/FlexDrivers/NewEVStatus/OptV6/' + str(month)
if not os.path.exists(directory):
os.makedirs(directory)
# run the daily optimization
for day in rng:
DailyPower = EVSE_Power[day]
DailyPower = DailyPower.replace('',np.nan)
DailyPower = DailyPower.dropna(axis=1,how='all')
DriverOnList = DailyPower.columns
N = len(DriverOnList)
opt_event = []
past_event = []
NonSlack_event = []
DailyPower['power'] = 0
DailyPower['NonControlledPower'] = DailyPower.fillna(0).sum(axis=1)
DailyPower['ts'] = range(96) # Add the tag of time step for each event id
if DailyPower.empty:
print 'EV Data on ' + day + ' is empty!'
DailyPower.to_csv(directory + '/' + day +'.csv')
else:
for eventid in DriverOnList[::-1]:
if (int(eventid) in EV_ID_Random) & (hasSlack(eventid,DailyPower)):
opt_event, LS_ref, RS_ref = flexEvent(eventid, DailyPower, past_event)
df2 = DailyPower[opt_event].fillna(0).sum(axis=1)
df2[LS_ref:96] = 0
DailyPower['power'] = 0
DailyPower['power'] = DailyPower[[x for x in past_event if x not in opt_event]].fillna(0).sum(axis=1) + df2
DailyPower = OptimizeCostV6(DailyPower, LS_ref, eventid, past_event, opt_event, onpeak, midpeak, peak)
print 'Optimization in ' + str(eventid)
else:
NonSlack_event.append(eventid)
opt_event, LS_ref, RS_ref = flexEvent(eventid, DailyPower, past_event)
df2 = DailyPower[opt_event].fillna(0).sum(axis=1)
df2[LS_ref:96] = 0
DailyPower['power'] = 0
DailyPower['power'] = DailyPower[[x for x in past_event if x not in opt_event]].fillna(0).sum(axis=1) + df2
DailyPower = OptimizeCostV6(DailyPower, LS_ref, eventid, past_event, opt_event, onpeak, midpeak, peak)
print str(eventid) + ' has no slack for charging load move!'
print 'Optimization in ' + str(eventid)
DailyPower['power'] = 0
DailyPower['power'] = DailyPower[past_event].fillna(0).sum(axis=1)
print 'Optimization on ' + day + ' is done successfully'
DailyPower.to_csv(directory +'/'+ day +'-V6.csv')
# EVSE_Power[StartDate:EndDate].to_csv('results/FlexDrivers/NewEVStatus/OptV5/'+'2015-'+str(month)+'Original.csv')
In [37]:
def OptimizeCostV1(DailyPower, opt_event, onpeak, midpeak, peak):
DayPower = DailyPower.fillna(0)
# Create matrix of slots
X = LpVariable.dicts("ChargingStatus", [(k, i, j) for k in range(len(opt_event)) for i in range(96) for j in range(96)], 0, 1, LpBinary)
costs=[]
for tim in DayPower.index:
costs.append(getCostValueForTime(tim.month,tim.strftime('%H:%M'),tim.isoweekday()))
cost={}
for k in xrange(len(opt_event)):
for i in xrange(96):
for j in xrange(96):
cost[k,i,j]=costs[j]
# Create model objective
model = pulp.LpProblem("SmartCharging",pulp.LpMinimize)
# Create model variables
min_curMidPeak = LpVariable("Min MidPeak",0)
min_curOnPeak = LpVariable("Min OnPeak",0)
min_curPeak = LpVariable("Min Peak",0)
tim = DayPower.index[0]
if (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [1,2,3,4,5]):
curMidPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(34,86)]
curPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33 + min_curMidPeak*0.13
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
for item in curMidPeak:
model += min_curMidPeak >= item
model += min_curMidPeak >= midpeak
elif (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [6,7]):
curPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [1,2,3,4,5]):
curMidPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in (range(34,48)+range(72,86))]
curPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
curOnPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(48,72)]
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curOnPeak*18.74 + min_curPeak*17.33 + min_curMidPeak*5.23
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
for item in curOnPeak:
model += min_curOnPeak >= item
model += min_curOnPeak >= onpeak
for item in curMidPeak:
model += min_curMidPeak >= item
model += min_curMidPeak >= midpeak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [6,7]):
curPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
else:
print "Error in month, day and time"
# Create model constrains
for k in xrange(len(opt_event)):
oldpower = DailyPower.loc[DailyPower[opt_event[k]] >= 0,[opt_event[k],'ts']]
newpower = DayPower.loc[DayPower[opt_event[k]] >= 0,[opt_event[k]]]
power2 = newpower[opt_event[k]].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
# Set the precedence and the session duration constraints
for i in xrange(LS):
model += (lpSum(X[k,i,j] for j in range(96)) == 0)
# Down zeros
for i in xrange(LS+M,96):
model += (lpSum(X[k,i,j] for j in range(96)) == 0)
# Left zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(LS)) == 0)
# Right zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(RS+1,96)) == 0)
# Add the constraints
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(LS,RS+1)) == 1)
for i in xrange(LS,LS+M-1):
# model += (lpSum(X[i,j]*(j+1) for j in range(LS,RS+1)) == lpSum(X[i+1,j]*(j+1) for j in range(LS,RS+1))-1)
model += (lpSum(X[k,i,j]*(j+1) for j in range(LS,RS+1)) <= lpSum(X[k,i+1,j]*(j+1) for j in range(LS,RS+1))-1) # discrete 15-min charging load blocks
model.solve(PULP_CBC_CMD(fracGap = 0.05, maxSeconds = 600, threads = None))
# model.solve(GLPK_CMD(options=['--mipgap','0.05']))
if model.status == 1:
for k in xrange(len(opt_event)):
oldpower = DailyPower.loc[DailyPower[opt_event[k]] >= 0,[opt_event[k],'ts']]
newpower = DayPower.loc[DayPower[opt_event[k]] >= 0,[opt_event[k]]]
power2 = newpower[opt_event[k]].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
slots = xrange(1,len(power2)+1)
sol=[]
for j in xrange(96):
sol.append(sum(X[k,i,j].varValue*(j) for i in xrange(LS,LS+M)))
corSlots=[slots[int(s)] for s in sol]
ct=0
power=[p*0 for p in power2]
for s in sol:
if s>0:
power[int(s)]=positive_power[ct]
ct=ct+1
else:
continue
DailyPower.loc[DailyPower.index[LS:RS+1],opt_event[k]] = power[LS:RS+1]
else:
print('The solution is suboptimal')
return DailyPower, model.status
In [39]:
# run the optimization for each month's EV charging sessions
noFeasibleDays = []
for month in range(7,8):
# monthly date range
StartDate = ts2015[ts2015['month'] == month].day[0]
EndDate = ts2015[ts2015['month'] == month].day[-1]
rng = map(lambda x: x.strftime('%Y-%m-%d'),pd.date_range(StartDate, EndDate, freq='D'))
EndDate2 = (datetime.strptime(EndDate,'%Y-%m-%d')+timedelta(days=1)).strftime('%Y-%m-%d')
# EVSE_Power = query_EVSE_eventids(StartDate,EndDate2)
# EVSE_Power = EVSE_Power[StartDate:EndDate]
EVSE_Power = pd.read_csv('results/FlexDrivers/NewEVStatus/2015EV/'+'2015-'+str(month)+'Original.csv')
EVSE_Power.index = EVSE_Power.ix[:,0].apply(lambda x: datetime.strptime(x,'%Y-%m-%d %H:%M:%S'))
EVSE_Power.drop(EVSE_Power.columns[[0]], axis=1, inplace=True)
# initialize the past peak demand values
DailyPower = []
peak=0
onpeak=0
midpeak=0
monthpeak = []
monthonpeak = []
monthmidpeak = []
directory = 'results/FlexDrivers/NewEVStatus/OptV1/' + str(month)
if not os.path.exists(directory):
os.makedirs(directory)
# run the daily optimization
for day in rng[0:1]:
DailyPower = EVSE_Power[day]
DailyPower = DailyPower.replace('',np.nan)
DailyPower = DailyPower.dropna(axis=1,how='all')
DriverOnList = DailyPower.columns
N = len(DriverOnList)
opt_event = []
DailyPower['power'] = 0
DailyPower['NonControlledPower'] = DailyPower.fillna(0).sum(axis=1)
DailyPower['ts'] = range(96) # Add the tag of time step for each event id
if DailyPower.empty:
print 'EV Data on ' + day + ' is empty!'
DailyPower.to_csv(directory + '/' + day +'.csv')
else:
for eventid in DriverOnList[::-1]:
# if (int(eventid) in EV_ID_Random) & (hasSlack(eventid,DailyPower)):
if (int(eventid) in EV_ID_Random):
opt_event.append(eventid)
DailyPower['power'] = DailyPower[[x for x in DriverOnList[::-1] if x not in opt_event]].fillna(0).sum(axis=1)
DailyPower,status = OptimizeCostV1(DailyPower, opt_event, onpeak, midpeak, peak)
if status == 1:
DailyPower['power'] = DailyPower['power'] + DailyPower[opt_event].fillna(0).sum(axis=1)
print 'Optimization on ' + day + ' is done successfully'
DailyPower.to_csv(directory +'/'+ day +'-V1.csv')
else:
copyfile('results/FlexDrivers/NewEVStatus/OptV3/' + str(month)+'/'+ day +'-V3.csv',directory +'/'+ day +'-V1.csv')
noFeasibleDays.append(day)
# DailyPower.to_csv(directory +'/'+ day +'-V1.csv')
# EVSE_Power[StartDate:EndDate].to_csv('results/FlexDrivers/NewEVStatus/OptV2/'+'2015-'+str(month)+'Original.csv')
In [17]:
with open(directory +'/'+'noFeasibleDays-V1.csv', 'wb') as f:
writer = csv.writer(f)
for val in noFeasibleDays:
writer.writerow([val])
In [18]:
def OptimizeCostV2(DailyPower, opt_event, onpeak, midpeak, peak):
DayPower = DailyPower.fillna(0)
# Create matrix of slots
X = LpVariable.dicts("ChargingStatus", [(k, i, j) for k in range(len(opt_event)) for i in range(96) for j in range(96)], 0, 1, LpBinary)
costs=[]
for tim in DayPower.index:
costs.append(getCostValueForTime(tim.month,tim.strftime('%H:%M'),tim.isoweekday()))
cost={}
for k in xrange(len(opt_event)):
for i in xrange(96):
for j in xrange(96):
cost[k,i,j]=costs[j]
# Create model objective
model = pulp.LpProblem("SmartCharging",pulp.LpMinimize)
# Create model variables
min_curMidPeak = LpVariable("Min MidPeak",0)
min_curOnPeak = LpVariable("Min OnPeak",0)
min_curPeak = LpVariable("Min Peak",0)
tim = DayPower.index[0]
if (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [1,2,3,4,5]):
curMidPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(34,86)]
curPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33 + min_curMidPeak*0.13
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
for item in curMidPeak:
model += min_curMidPeak >= item
model += min_curMidPeak >= midpeak
elif (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [6,7]):
curPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [1,2,3,4,5]):
curMidPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in (range(34,48)+range(72,86))]
curPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
curOnPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(48,72)]
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curOnPeak*18.74 + min_curPeak*17.33 + min_curMidPeak*5.23
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
for item in curOnPeak:
model += min_curOnPeak >= item
model += min_curOnPeak >= onpeak
for item in curMidPeak:
model += min_curMidPeak >= item
model += min_curMidPeak >= midpeak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [6,7]):
curPeak = [lpSum(lpSum(DayPower[opt_event[k]].tolist()[i]*X[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[opt_event[k]].tolist()[i])*X[k,i,j]*cost[k,i,j] for k in xrange(len(opt_event)) for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
else:
print "Error in month, day and time"
# Create model constrains
for k in xrange(len(opt_event)):
oldpower = DailyPower.loc[DailyPower[opt_event[k]] >= 0,[opt_event[k],'ts']]
newpower = DayPower.loc[DayPower[opt_event[k]] >= 0,[opt_event[k]]]
power2 = newpower[opt_event[k]].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
# Set the precedence and the session duration constraints
for i in xrange(LS):
model += (lpSum(X[k,i,j] for j in range(96)) == 0)
# Down zeros
for i in xrange(LS+M,96):
model += (lpSum(X[k,i,j] for j in range(96)) == 0)
# Left zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(LS)) == 0)
# Right zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(RS+1,96)) == 0)
# Add the constraints
for i in xrange(LS,LS+M):
model += (lpSum(X[k,i,j] for j in range(LS,RS+1)) == 1)
for i in xrange(LS,LS+M-1):
model += (lpSum(X[i,j]*(j+1) for j in range(LS,RS+1)) == lpSum(X[i+1,j]*(j+1) for j in range(LS,RS+1))-1)
# model += (lpSum(X[k,i,j]*(j+1) for j in range(LS,RS+1)) <= lpSum(X[k,i+1,j]*(j+1) for j in range(LS,RS+1))-1) # discrete 15-min charging load blocks
model.solve(PULP_CBC_CMD(fracGap = 0.05, maxSeconds = 100, threads = None))
# model.solve(GLPK_CMD(options=['--mipgap','0.05']))
if model.status == 1:
for k in xrange(len(opt_event)):
oldpower = DailyPower.loc[DailyPower[opt_event[k]] >= 0,[opt_event[k],'ts']]
newpower = DayPower.loc[DayPower[opt_event[k]] >= 0,[opt_event[k]]]
power2 = newpower[opt_event[k]].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
slots = xrange(1,len(power2)+1)
sol=[]
for j in xrange(96):
sol.append(sum(X[k,i,j].varValue*(j) for i in xrange(LS,LS+M)))
corSlots=[slots[int(s)] for s in sol]
ct=0
power=[p*0 for p in power2]
for s in sol:
if s>0:
power[int(s)]=positive_power[ct]
ct=ct+1
else:
continue
DailyPower.loc[DailyPower.index[LS:RS+1],opt_event[k]] = power[LS:RS+1]
else:
print('The solution is suboptimal')
return DailyPower, model.status
In [19]:
# run the optimization for each month's EV charging sessions
noFeasibleDays = []
for month in range(1,13):
# monthly date range
StartDate = ts2015[ts2015['month'] == month].day[0]
EndDate = ts2015[ts2015['month'] == month].day[-1]
rng = map(lambda x: x.strftime('%Y-%m-%d'),pd.date_range(StartDate, EndDate, freq='D'))
EndDate2 = (datetime.strptime(EndDate,'%Y-%m-%d')+timedelta(days=1)).strftime('%Y-%m-%d')
# EVSE_Power = query_EVSE_eventids(StartDate,EndDate2)
# EVSE_Power = EVSE_Power[StartDate:EndDate]
EVSE_Power = pd.read_csv('results/FlexDrivers/NewEVStatus/2015EV/'+'2015-'+str(month)+'Original.csv')
EVSE_Power.index = EVSE_Power.ix[:,0].apply(lambda x: datetime.strptime(x,'%Y-%m-%d %H:%M:%S'))
EVSE_Power.drop(EVSE_Power.columns[[0]], axis=1, inplace=True)
# initialize the past peak demand values
DailyPower = []
peak=0
onpeak=0
midpeak=0
monthpeak = []
monthonpeak = []
monthmidpeak = []
directory = 'results/FlexDrivers/NewEVStatus/OptV2/' + str(month)
if not os.path.exists(directory):
os.makedirs(directory)
# run the daily optimization
for day in rng:
DailyPower = EVSE_Power[day]
DailyPower = DailyPower.replace('',np.nan)
DailyPower = DailyPower.dropna(axis=1,how='all')
DriverOnList = DailyPower.columns
N = len(DriverOnList)
opt_event = []
DailyPower['power'] = 0
DailyPower['NonControlledPower'] = DailyPower.fillna(0).sum(axis=1)
DailyPower['ts'] = range(96) # Add the tag of time step for each event id
if DailyPower.empty:
print 'EV Data on ' + day + ' is empty!'
DailyPower.to_csv(directory + '/' + day +'.csv')
else:
for eventid in DriverOnList[::-1]:
if (int(eventid) in EV_ID_Random) & (hasSlack(eventid,DailyPower)):
opt_event.append(eventid)
DailyPower['power'] = DailyPower[[x for x in DriverOnList[::-1] if x not in opt_event]].fillna(0).sum(axis=1)
DailyPower,status = OptimizeCostV2(DailyPower, opt_event, onpeak, midpeak, peak)
if status == 1:
DailyPower['power'] = DailyPower['power'] + DailyPower[opt_event].fillna(0).sum(axis=1)
print 'Optimization on ' + day + ' is done successfully'
DailyPower.to_csv(directory +'/'+ day +'-V2.csv')
else:
copyfile('results/FlexDrivers/NewEVStatus/OptV3/' + str(month)+'/'+ day +'-V3.csv',directory +'/'+ day +'-V2.csv')
noFeasibleDays.append(day)
# DailyPower.to_csv(directory +'/'+ day +'-V1.csv')
# EVSE_Power[StartDate:EndDate].to_csv('results/FlexDrivers/NewEVStatus/OptV2/'+'2015-'+str(month)+'Original.csv')
with open(directory +'/'+'noFeasibleDays-V2.csv', 'wb') as f:
writer = csv.writer(f)
for val in noFeasibleDays:
writer.writerow([val])
In [10]:
def OptimizeCostV3(DailyPower, eventid, onpeak, midpeak, peak):
# DailyPower2 = DailyPower.copy()
# DailyPower2.loc[DailyPower2.index[0:LS_ref],opt_event] = np.nan
DayPower = DailyPower.fillna(0)
# Create matrix of slots
X = LpVariable.dicts("ChargingStatus", [(i, j) for i in range(96) for j in range(96)], 0, 1, LpBinary)
costs=[]
for tim in DayPower.index:
costs.append(getCostValueForTime(tim.month,tim.strftime('%H:%M'),tim.isoweekday()))
cost={}
for i in xrange(96):
for j in xrange(96):
cost[i,j]=costs[j]
# Create model objective
model = pulp.LpProblem("SmartCharging",pulp.LpMinimize)
# Create model variables
min_curMidPeak = LpVariable("Min MidPeak",0)
min_curOnPeak = LpVariable("Min OnPeak",0)
min_curPeak = LpVariable("Min Peak",0)
tim = DayPower.index[0]
if (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [1,2,3,4,5]):
curMidPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(34,86)]
curPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[eventid].tolist()[i])*X[i,j]*cost[i,j] for i in xrange(96) for j in xrange(96)) + min_curPeak*17.33 + min_curMidPeak*0.13
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
for item in curMidPeak:
model += min_curMidPeak >= item
model += min_curMidPeak >= midpeak
elif (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [6,7]):
curPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[eventid].tolist()[i])*X[i,j]*cost[i,j] for i in xrange(96) for j in xrange(96)) + min_curPeak*17.33
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [1,2,3,4,5]):
curMidPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in (range(34,48)+range(72,86))]
curPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
curOnPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(48,72)]
model += lpSum((DayPower[eventid].tolist()[i])*X[i,j]*cost[i,j] for i in xrange(96) for j in xrange(96)) + \
min_curOnPeak*18.74 + min_curPeak*17.33 + min_curMidPeak*5.23
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
for item in curOnPeak:
model += min_curOnPeak >= item
model += min_curOnPeak >= onpeak
for item in curMidPeak:
model += min_curMidPeak >= item
model += min_curMidPeak >= midpeak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [6,7]):
curPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[eventid].tolist()[i])*X[i,j]*cost[i,j] for i in xrange(96) for j in xrange(96)) + min_curPeak*17.33
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
else:
print "Error in month, day and time"
# Create model constrains
oldpower = DailyPower.loc[DailyPower[eventid] >= 0,[eventid,'ts']]
newpower = DayPower.loc[DayPower[eventid] >= 0,[eventid]]
power2 = newpower[eventid].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
# Set the precedence and the session duration constraints
for i in xrange(LS):
model += (lpSum(X[i,j] for j in range(96)) == 0)
# Down zeros
for i in xrange(LS+M,96):
model += (lpSum(X[i,j] for j in range(96)) == 0)
# Left zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[i,j] for j in range(LS)) == 0)
# Right zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[i,j] for j in range(RS+1,96)) == 0)
# Add the constraints
for i in xrange(LS,LS+M):
model += (lpSum(X[i,j] for j in range(LS,RS+1)) == 1)
for i in xrange(LS,LS+M-1):
# model += (lpSum(X[i,j]*(j+1) for j in range(LS,RS+1)) == lpSum(X[i+1,j]*(j+1) for j in range(LS,RS+1))-1)
model += (lpSum(X[i,j]*(j+1) for j in range(LS,RS+1)) <= lpSum(X[i+1,j]*(j+1) for j in range(LS,RS+1))-1) # discrete 15-min charging load blocks
model.solve()
if model.status == 1:
oldpower = DailyPower.loc[DailyPower[eventid] >= 0,[eventid,'ts']]
newpower = DayPower.loc[DayPower[eventid] >= 0,[eventid]]
power2 = newpower[eventid].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
slots = xrange(1,len(power2)+1)
sol=[]
for j in xrange(96):
sol.append(sum(X[i,j].varValue*(j) for i in xrange(LS,LS+M)))
corSlots=[slots[int(s)] for s in sol]
ct=0
power=[p*0 for p in power2]
for s in sol:
if s>0:
power[int(s)]=positive_power[ct]
ct=ct+1
else:
continue
DailyPower.loc[DailyPower.index[LS:RS+1],eventid] = power[LS:RS+1]
else:
print('The solution is suboptimal')
return DailyPower
In [11]:
# run the optimization for each month's EV charging sessions
for month in range(3,13):
# monthly date range
StartDate = ts2015[ts2015['month'] == month].day[0]
EndDate = ts2015[ts2015['month'] == month].day[-1]
rng = map(lambda x: x.strftime('%Y-%m-%d'),pd.date_range(StartDate, EndDate, freq='D'))
EndDate2 = (datetime.strptime(EndDate,'%Y-%m-%d')+timedelta(days=1)).strftime('%Y-%m-%d')
# EVSE_Power = query_EVSE_eventids(StartDate,EndDate2)
# EVSE_Power = EVSE_Power[StartDate:EndDate]
EVSE_Power = pd.read_csv('results/FlexDrivers/NewEVStatus/2015EV/'+'2015-'+str(month)+'Original.csv')
EVSE_Power.index = EVSE_Power.ix[:,0].apply(lambda x: datetime.strptime(x,'%Y-%m-%d %H:%M:%S'))
EVSE_Power.drop(EVSE_Power.columns[[0]], axis=1, inplace=True)
# initialize the past peak demand values
DailyPower = []
peak=0
onpeak=0
midpeak=0
monthpeak = []
monthonpeak = []
monthmidpeak = []
directory = 'results/FlexDrivers/NewEVStatus/OptV3/' + str(month)
if not os.path.exists(directory):
os.makedirs(directory)
# run the daily optimization
for day in rng:
DailyPower = EVSE_Power[day]
DailyPower = DailyPower.replace('',np.nan)
DailyPower = DailyPower.dropna(axis=1,how='all')
DriverOnList = DailyPower.columns
N = len(DriverOnList)
opt_event = []
past_event = []
DailyPower['power'] = 0
DailyPower['NonControlledPower'] = DailyPower[DriverOnList].fillna(0).sum(axis=1)
DailyPower['ts'] = range(96) # Add the tag of time step for each event id
if DailyPower.empty:
print 'EV Data on ' + day + ' is empty!'
DailyPower.to_csv(directory + '/' + day +'.csv')
else:
for eventid in DriverOnList[::-1]:
if (int(eventid) in EV_ID_Random) & (hasSlack(eventid,DailyPower)):
# Run the optimization for the current EV
DailyPower = OptimizeCostV3(DailyPower, eventid, onpeak, midpeak, peak)
else:
print 'No slack in '+str(eventid)
DailyPower['power'] = DailyPower['power'] + DailyPower[eventid].fillna(0)
print 'Optimization on ' + day + ' is done successfully'
DailyPower.to_csv(directory +'/'+ day +'-V3.csv')
# EVSE_Power[StartDate:EndDate].to_csv('results/FlexDrivers/NewEVStatus/OptV4/'+'2015-'+str(month)+'Original.csv')
In [ ]:
def OptimizeCostV4(DailyPower, eventid, onpeak, midpeak, peak):
# DailyPower2 = DailyPower.copy()
# DailyPower2.loc[DailyPower2.index[0:LS_ref],opt_event] = np.nan
DayPower = DailyPower.fillna(0)
# Create matrix of slots
X = LpVariable.dicts("ChargingStatus", [(i, j) for i in range(96) for j in range(96)], 0, 1, LpBinary)
costs=[]
for tim in DayPower.index:
costs.append(getCostValueForTime(tim.month,tim.strftime('%H:%M'),tim.isoweekday()))
cost={}
for i in xrange(96):
for j in xrange(96):
cost[i,j]=costs[j]
# Create model objective
model = pulp.LpProblem("SmartCharging",pulp.LpMinimize)
# Create model variables
min_curMidPeak = LpVariable("Min MidPeak",0)
min_curOnPeak = LpVariable("Min OnPeak",0)
min_curPeak = LpVariable("Min Peak",0)
tim = DayPower.index[0]
if (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [1,2,3,4,5]):
curMidPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(34,86)]
curPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[eventid].tolist()[i])*X[i,j]*cost[i,j] for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33 + min_curMidPeak*0.13
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
for item in curMidPeak:
model += min_curMidPeak >= item
model += min_curMidPeak >= midpeak
elif (tim.month in [1,2,3,4,11,12]) & (tim.isoweekday() in [6,7]):
curPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[eventid].tolist()[i])*X[i,j]*cost[i,j] for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [1,2,3,4,5]):
curMidPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in (range(34,48)+range(72,86))]
curPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
curOnPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(48,72)]
model += lpSum((DayPower[eventid].tolist()[i])*X[i,j]*cost[i,j] for i in xrange(96) for j in xrange(96)) + \
min_curOnPeak*18.74 + min_curPeak*17.33 + min_curMidPeak*5.23
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
for item in curOnPeak:
model += min_curOnPeak >= item
model += min_curOnPeak >= onpeak
for item in curMidPeak:
model += min_curMidPeak >= item
model += min_curMidPeak >= midpeak
elif (tim.month in [5,6,7,8,9,10]) & (tim.isoweekday() in [6,7]):
curPeak = [lpSum(lpSum(DayPower[eventid].tolist()[i]*X[i,j] for i in xrange(96)) + DayPower['power'].tolist()[j]) for j in xrange(96)]
model += lpSum((DayPower[eventid].tolist()[i])*X[i,j]*cost[i,j] for i in xrange(96) for j in xrange(96)) + \
min_curPeak*17.33
for item in curPeak:
model += min_curPeak >= item
model += min_curPeak >= peak
else:
print "Error in month, day and time"
# Create model constrains
oldpower = DailyPower.loc[DailyPower[eventid] >= 0,[eventid,'ts']]
newpower = DayPower.loc[DayPower[eventid] >= 0,[eventid]]
power2 = newpower[eventid].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
# Set the precedence and the session duration constraints
for i in xrange(LS):
model += (lpSum(X[i,j] for j in range(96)) == 0)
# Down zeros
for i in xrange(LS+M,96):
model += (lpSum(X[i,j] for j in range(96)) == 0)
# Left zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[i,j] for j in range(LS)) == 0)
# Right zeros
for i in xrange(LS,LS+M):
model += (lpSum(X[i,j] for j in range(RS+1,96)) == 0)
# Add the constraints
for i in xrange(LS,LS+M):
model += (lpSum(X[i,j] for j in range(LS,RS+1)) == 1)
for i in xrange(LS,LS+M-1):
model += (lpSum(X[i,j]*(j+1) for j in range(LS,RS+1)) == lpSum(X[i+1,j]*(j+1) for j in range(LS,RS+1))-1)
# model += (lpSum(X[k,i,j]*(j+1) for j in range(LS,RS+1)) <= lpSum(X[k,i+1,j]*(j+1) for j in range(LS,RS+1))-1) # discrete 15-min charging load blocks
model.solve()
if model.status == 1:
oldpower = DailyPower.loc[DailyPower[eventid] >= 0,[eventid,'ts']]
newpower = DayPower.loc[DayPower[eventid] >= 0,[eventid]]
power2 = newpower[eventid].tolist()
positive_power = [p for p in power2 if p>0]
M=len(positive_power)
# Left and Right slots
LS = oldpower['ts'][0]
RS = oldpower['ts'][-1]
slots = xrange(1,len(power2)+1)
sol=[]
for j in xrange(96):
sol.append(sum(X[i,j].varValue*(j) for i in xrange(LS,LS+M)))
corSlots=[slots[int(s)] for s in sol]
ct=0
power=[p*0 for p in power2]
for s in sol:
if s>0:
power[int(s)]=positive_power[ct]
ct=ct+1
else:
continue
DailyPower.loc[DailyPower.index[LS:RS+1],eventid] = power[LS:RS+1]
else:
print('The solution is suboptimal')
return DailyPower
In [ ]:
# run the optimization for each month's EV charging sessions
for month in range(1,13):
# monthly date range
StartDate = ts2015[ts2015['month'] == month].day[0]
EndDate = ts2015[ts2015['month'] == month].day[-1]
rng = map(lambda x: x.strftime('%Y-%m-%d'),pd.date_range(StartDate, EndDate, freq='D'))
EndDate2 = (datetime.strptime(EndDate,'%Y-%m-%d')+timedelta(days=1)).strftime('%Y-%m-%d')
# EVSE_Power = query_EVSE_eventids(StartDate,EndDate2)
# EVSE_Power = EVSE_Power[StartDate:EndDate]
EVSE_Power = pd.read_csv('results/FlexDrivers/NewEVStatus/2015EV/'+'2015-'+str(month)+'Original.csv')
EVSE_Power.index = EVSE_Power.ix[:,0].apply(lambda x: datetime.strptime(x,'%Y-%m-%d %H:%M:%S'))
EVSE_Power.drop(EVSE_Power.columns[[0]], axis=1, inplace=True)
# initialize the past peak demand values
DailyPower = []
peak=0
onpeak=0
midpeak=0
monthpeak = []
monthonpeak = []
monthmidpeak = []
directory = 'results/FlexDrivers/NewEVStatus/OptV4/' + str(month)
if not os.path.exists(directory):
os.makedirs(directory)
# run the daily optimization
for day in rng:
DailyPower = EVSE_Power[day]
DailyPower = DailyPower.replace('',np.nan)
DailyPower = DailyPower.dropna(axis=1,how='all')
DriverOnList = DailyPower.columns
N = len(DriverOnList)
opt_event = []
past_event = []
DailyPower['power'] = 0
DailyPower['NonControlledPower'] = DailyPower[DriverOnList].fillna(0).sum(axis=1)
DailyPower['ts'] = range(96) # Add the tag of time step for each event id
if DailyPower.empty:
print 'EV Data on ' + day + ' is empty!'
DailyPower.to_csv(directory + '/' + day +'.csv')
else:
for eventid in DriverOnList[::-1]:
if (int(eventid) in EV_ID_Random) & (hasSlack(eventid,DailyPower)):
# Run the optimization for the current EV
DailyPower = OptimizeCostV4(DailyPower, eventid, onpeak, midpeak, peak)
else:
print 'No slack in '+str(eventid)
DailyPower['power'] = DailyPower['power'] + DailyPower[eventid].fillna(0)
print 'Optimization on ' + day + ' is done successfully'
DailyPower.to_csv(directory +'/'+ day +'-V4.csv')
# EVSE_Power[StartDate:EndDate].to_csv('results/FlexDrivers/NewEVStatus/OptV4/'+'2015-'+str(month)+'Original.csv')
In [ ]: