In [12]:
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
import timeit


from IPython.display import clear_output
from IPython.display import HTML
from matplotlib.ticker import MultipleLocator

Meter data at ALCOPARK

Read the utility rate structure of E19, PGE


In [13]:
# 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 = utility['Secondary'][1]
summer_midpeak_demand = utility['Secondary'][2]
summer_monthly_demand = utility['Secondary'][3]
winter_midpeak_demand = utility['Secondary'][4]
winter_monthly_demand = utility['Secondary'][5]

# energy charge
summer_onpeak_energy = utility['Secondary'][6]
summer_midpeak_energy = utility['Secondary'][7]
summer_offpeak_energy = utility['Secondary'][8]
winter_midpeak_energy = utility['Secondary'][9]
winter_offpeak_energy = 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]

Read original Nonfleet EVs and Optimized EVs


In [84]:
alcopark = pd.read_csv('results/FlexDrivers/NewEVStatus/alcopark/' + '2015-EVPower.csv')
alcopark.index = alcopark.ix[:,0].apply(lambda x: datetime.strptime(x,'%m/%d/%y %H:%M'))
alcopark.drop(alcopark.columns[[0]], axis=1, inplace=True)
alcopark.index.name = 'time'

alcopark['day'] = map(lambda x: x.strftime('%Y-%m-%d'),alcopark.index)
alcopark['weekday'] = alcopark.index.weekday
alcopark['hour'] = alcopark.index.hour
alcopark['Date']=alcopark.index
alcopark['time'] = alcopark.Date.apply(lambda x: x.strftime('%H:%M'))
alcopark['month'] = alcopark.index.month

alcopark['basepower'] = alcopark.building - alcopark.power
alcopark['hour'] = alcopark.hour.astype(int)
alcopark['month'] = alcopark.month.astype(int)
alcopark['holiday'] = alcopark.Date.apply(lambda x: x.strftime('%Y-%m-%d') in holidays.day.values)
alcopark['holiday'] = alcopark.holiday.astype(int)
alcopark['PeakPeriod'] = 0

alcopark.loc[alcopark['holiday']==1,'weekday'] = -1

alcopark.loc[(alcopark.month >= 5) & (alcopark.month <= 10) & \
             (alcopark.weekday >= 0) & (alcopark.weekday <= 4) & \
             (alcopark.time >= '08:30') & (alcopark.time < '12:00'),'PeakPeriod'] = 1
alcopark.loc[(alcopark.month >= 5) & (alcopark.month <= 10) & \
             (alcopark.weekday >= 0) & (alcopark.weekday <= 4) & \
             (alcopark.time >= '12:00') & (alcopark.time < '18:00'),'PeakPeriod'] = 2
alcopark.loc[(alcopark.month >= 5) & (alcopark.month <= 10) & \
             (alcopark.weekday >= 0) & (alcopark.weekday <= 4) & \
             (alcopark.time >= '18:00') & (alcopark.time < '21:30'),'PeakPeriod'] = 3

alcopark.loc[(alcopark.month.apply(lambda x: (x in [1,2,3,4,11,12]))) & \
             (alcopark.weekday >= 0) & (alcopark.weekday <= 4) & \
             (alcopark.time >= '08:30') & (alcopark.time < '21:30'),'PeakPeriod'] = 1

In [45]:
# Plot Arrival and departure time characteristics.
sns.set_context('poster')
fig, (ax1, ax2) = plt.subplots(1,2, figsize=(18, 6), facecolor='w', edsgecolor='k', sharex=False)
# fig.subplots_adjust(hspace = .1, wspace=.001)

num_bins = range(24)
ax1.scatter(Alco['2015-07-29'].hour, Alco['2015-07-29'].building,alpha=0.5)
ax2.scatter(Alco['2015-07-29'].hour, Alco['2015-07-29'].basepower, facecolor='red', alpha=0.5)
ax1.set_ylabel('Power, [kW]')
ax1.set_xlabel('Hour of Day')
ax1.set_title('Whole Building Power at Alcopark')
ax1.set_xlim(0,24)

ax2.set_ylabel('Power, [kW]')
ax2.set_xlabel('Hour of Day')
ax2.set_title('Base Power at Alcopark')
ax2.set_xlim(0,24)


Out[45]:
(0, 24)

In [ ]:
def alcopark_power(alcopark, month):
    df_wd = alcopark.loc[(alcopark['month']==month) & (alcopark['weekday'] >= 0) & (alcopark['weekday'] <= 4)]
    df_wk = alcopark.loc[(alcopark['month']==month) & (alcopark.weekday.apply(lambda x: (x in [5,6])))]

    df_bldg_wd = df_wd.pivot(index='time',columns='day',values='building')
    df_base_wd = df_wd.pivot(index='time',columns='day',values='basepower')
    df_ev_wd = df_wd.pivot(index='time',columns='day',values='power')

    df_bldg_wk = df_wk.pivot(index='time',columns='day',values='building')
    df_base_wk = df_wk.pivot(index='time',columns='day',values='basepower')

    fig, (axes) = plt.subplots(2,2, figsize=(18, 12), facecolor='w', edgecolor='k', sharey=True, sharex=True)
    fig.subplots_adjust(hspace = .1, wspace=.1)

    df_bldg_wd.plot(ax=axes[0,0],c='gray',alpha=0.8)
    axes[0,0].plot(df_bldg_wd.mean(axis=1),c='b',label='Average Meter Power')
    axes[0,0].set_ylabel('Power, [kW]')
    axes[0,0].set_xlabel('Hour of Day')
    axes[0,0].set_title('Whole Building Power on Weekdays in '+'2015/'+str(month))
    axes[0,0].legend().set_visible(False)

    df_base_wd.plot(ax=axes[0,1],c='gray',alpha=0.8)
    axes[0,1].plot(df_base_wd.mean(axis=1),c='red', label='Average Base Power')
    axes[0,1].set_ylabel('Power, [kW]')
    axes[0,1].set_xlabel('Hour of Day')
    axes[0,1].set_title('Base Power without EVs on Weekdays in '+'2015/'+str(month))
    axes[0,1].legend().set_visible(False)

    df_ev_wd.plot(ax=axes[1,0],c='gray',alpha=0.8)
    axes[1,0].plot(df_ev_wd.mean(axis=1),color='g')
    axes[1,0].set_ylabel('Power, [kW]')
    axes[1,0].set_xlabel('Hour of Day')
    axes[1,0].set_title('EVs Power on Weekdays in '+'2015/'+str(month))
    axes[1,0].legend().set_visible(False)

    # df_base_wk.plot(ax=axes[1,1],color='gray',alpha=0.8)
    axes[1,1].plot(df_bldg_wd.mean(axis=1),c='b',label='Meter Power')
    axes[1,1].plot(df_base_wd.mean(axis=1),c='r',label='Base Power')

    axes[1,1].set_ylabel('Power, [kW]')
    axes[1,1].set_xlabel('Hour of Day')
    axes[1,1].set_title('Base Power, EVs on Weekdays in '+'2015/'+str(month))
    axes[1,1].legend()
    axes[1,1].fill_between(range(len(df_base_wd.mean(axis=1))), df_base_wd.mean(axis=1), df_bldg_wd.mean(axis=1), facecolor='gray')

    fig.savefig('results/FlexDrivers/NewEVStatus/alcopark/Alcopark-'+str(month)+'.png', dpi=300, format='png',bbox_inches='tight')

In [53]:
sns.set_context('poster')
for month in range(1,13):
    alcopark_power(alcopark,month)
    print 'month ' + str(month) + ' is done.'


month 1 is done.
month 2 is done.
month 3 is done.
month 4 is done.
month 5 is done.
month 6 is done.
month 7 is done.
month 8 is done.
month 9 is done.
month 10 is done.
month 11 is done.
month 12 is done.

Peak load, load variability


In [29]:
month=1
df_wd = alcopark.loc[(alcopark['month']==month) & (alcopark['weekday'] >= 0) & (alcopark['weekday'] <= 4)]
df_wk = alcopark.loc[(alcopark['month']==month) & (alcopark.weekday.apply(lambda x: (x in [5,6])))]

df_bldg_wd = df_wd.pivot(index='time',columns='day',values='building')
df_base_wd = df_wd.pivot(index='time',columns='day',values='basepower')
df_ev_wd = df_wd.pivot(index='time',columns='day',values='power')

df_bldg_wk = df_wk.pivot(index='time',columns='day',values='building')
df_base_wk = df_wk.pivot(index='time',columns='day',values='basepower')

In [83]:
fig, (ax1, ax2) = plt.subplots(1,2, figsize=(18, 6), facecolor='w', edgecolor='k', sharey=True)
fig.subplots_adjust(hspace = .1, wspace=.1)

# num_bins = range(24)
ax1.plot(df_bldg_wd.std(axis=1)/df_bldg_wd.mean(axis=1))
ax2.plot(df_base_wd.std(axis=1)/df_base_wd.mean(axis=1))
ax1.set_ylabel('Power, [kW]')
ax1.set_xlabel('Hour of Day')
ax1.set_title('Whole Building Power Variability at Alcopark')
ax1.set_xticklabels(range(0,96,4))

majorLocator   = MultipleLocator(16)
ax2.xaxis.set_major_locator(majorLocator)
ax2.set_ylabel('Power, [kW]')
ax2.set_xlabel('Hour of Day')
ax2.set_title('Base Power Variability at Alcopark')
ax2.set_xlim(0,96)
ax2.set_xticklabels(range(-4,96,4))


Out[83]:
[<matplotlib.text.Text at 0x115e2ead0>,
 <matplotlib.text.Text at 0x13d0015d0>,
 <matplotlib.text.Text at 0x138f7cf10>,
 <matplotlib.text.Text at 0x138f72a10>,
 <matplotlib.text.Text at 0x138f79990>,
 <matplotlib.text.Text at 0x138f702d0>,
 <matplotlib.text.Text at 0x138f49650>,
 <matplotlib.text.Text at 0x11502c990>,
 <matplotlib.text.Text at 0x11502cf90>]

In [117]:
fig, axes = plt.subplots(1,3, figsize=(21, 5), facecolor='w', edgecolor='k', sharey=True)
fig.subplots_adjust(hspace = .2, wspace=.1)

for mon in [5,6,7]:
    df_wd = alcopark.loc[(alcopark['month']==mon)]
    df_peaks = df_wd.pivot(index='Date',columns='PeakPeriod',values='building')

    df_peaks.boxplot(ax=axes[mon-5])
    axes[mon-5].set_ylabel('Power, [kW]')
    axes[mon-5].set_title('Peak Power Distribution in 2015/'+str(mon))
    axes[mon-5].set_xticklabels(['Off-Peak','Mid-Peak','On-Peak','Mid-Peak'])
    axes[mon-5].set_ylim(0,120)

fig.savefig('results/FlexDrivers/NewEVStatus/alcopark/Alcopark-Peak-1.png', dpi=300, format='png',bbox_inches='tight')

fig, axes = plt.subplots(1,3, figsize=(21, 5), facecolor='w', edgecolor='k', sharey=True)
fig.subplots_adjust(hspace = .2, wspace=.1)

for mon in [8,9,10]:
    df_wd = alcopark.loc[(alcopark['month']==mon)]
    df_peaks = df_wd.pivot(index='Date',columns='PeakPeriod',values='building')

    df_peaks.boxplot(ax=axes[mon-8])
    axes[mon-8].set_ylabel('Power, [kW]')
    axes[mon-8].set_title('Peak Power Distribution in 2015/'+str(mon))
    axes[mon-8].set_xticklabels(['Off-Peak','Mid-Peak','On-Peak','Mid-Peak'])
    axes[mon-8].set_ylim(0,120)

fig.savefig('results/FlexDrivers/NewEVStatus/alcopark/Alcopark-Peak-2.png', dpi=300, format='png',bbox_inches='tight')



In [118]:
fig, axes = plt.subplots(1,3, figsize=(21, 5), facecolor='w', edgecolor='k', sharey=True)
fig.subplots_adjust(hspace = .2, wspace=.1)

for mon in [5,6,7]:
    df_wd = alcopark.loc[(alcopark['month']==mon)]
    df_peaks = df_wd.pivot(index='Date',columns='PeakPeriod',values='power')

    df_peaks.boxplot(ax=axes[mon-5])
    axes[mon-5].set_ylabel('Power, [kW]')
    axes[mon-5].set_title('EV Peak Power Distribution in 2015/'+str(mon))
    axes[mon-5].set_xticklabels(['Off-Peak','Mid-Peak','On-Peak','Mid-Peak'])
#     axes[mon-5].set_ylim(0,120)

fig.savefig('results/FlexDrivers/NewEVStatus/alcopark/Alcopark-EVPeak-1.png', dpi=300, format='png',bbox_inches='tight')

fig, axes = plt.subplots(1,3, figsize=(21, 5), facecolor='w', edgecolor='k', sharey=True)
fig.subplots_adjust(hspace = .2, wspace=.1)

for mon in [8,9,10]:
    df_wd = alcopark.loc[(alcopark['month']==mon)]
    df_peaks = df_wd.pivot(index='Date',columns='PeakPeriod',values='power')

    df_peaks.boxplot(ax=axes[mon-8])
    axes[mon-8].set_ylabel('Power, [kW]')
    axes[mon-8].set_title('EV Peak Power Distribution in 2015/'+str(mon))
    axes[mon-8].set_xticklabels(['Off-Peak','Mid-Peak','On-Peak','Mid-Peak'])
#     axes[mon-8].set_ylim(0,120)

fig.savefig('results/FlexDrivers/NewEVStatus/alcopark/Alcopark-EVPeak-2.png', dpi=300, format='png',bbox_inches='tight')



In [119]:
fig, axes = plt.subplots(1,3, figsize=(21, 5), facecolor='w', edgecolor='k', sharey=True)
fig.subplots_adjust(hspace = .2, wspace=.1)

for mon in [5,6,7]:
    df_wd = alcopark.loc[(alcopark['month']==mon)]
    df_peaks = df_wd.pivot(index='Date',columns='PeakPeriod',values='basepower')

    df_peaks.boxplot(ax=axes[mon-5])
    axes[mon-5].set_ylabel('Power, [kW]')
    axes[mon-5].set_title('EV Peak Power Distribution in 2015/'+str(mon))
    axes[mon-5].set_xticklabels(['Off-Peak','Mid-Peak','On-Peak','Mid-Peak'])
#     axes[mon-5].set_ylim(0,120)

fig.savefig('results/FlexDrivers/NewEVStatus/alcopark/Alcopark-BasePeak-1.png', dpi=300, format='png',bbox_inches='tight')

fig, axes = plt.subplots(1,3, figsize=(21, 5), facecolor='w', edgecolor='k', sharey=True)
fig.subplots_adjust(hspace = .2, wspace=.1)

for mon in [8,9,10]:
    df_wd = alcopark.loc[(alcopark['month']==mon)]
    df_peaks = df_wd.pivot(index='Date',columns='PeakPeriod',values='basepower')

    df_peaks.boxplot(ax=axes[mon-8])
    axes[mon-8].set_ylabel('Power, [kW]')
    axes[mon-8].set_title('EV Peak Power Distribution in 2015/'+str(mon))
    axes[mon-8].set_xticklabels(['Off-Peak','Mid-Peak','On-Peak','Mid-Peak'])
#     axes[mon-8].set_ylim(0,120)

fig.savefig('results/FlexDrivers/NewEVStatus/alcopark/Alcopark-BasePeak-2.png', dpi=300, format='png',bbox_inches='tight')



In [ ]: