In [2]:
import pandas as pd
import numpy as np
import datetime
import plotly
import plotly.offline as po
import plotly.plotly as py
import plotly.graph_objs as go
from plotly.offline import download_plotlyjs
import cufflinks as cf
import cython
from netCDF4 import *
cf.go_offline()
cf.set_config_file(theme = 'pearl')
#Obtain free Plotly account to authenticate with plotly before entering off-line mode
#plotly.tools.set_credentials_file(username='PLOTLY_USERNAME', api_key='PLOTLY_API_KEY')
np.set_printoptions(suppress = True, precision = 3, linewidth = 100)
pd.options.display.float_format = '{:.3f}'.format
%load_ext cython
print("---Imports & Plotly Connection Complete---")
In [3]:
days = 365 # Days input for Leap Year handling
date_start = '01/01/2013'
date_end = '12/31/2013'
simulation_years = 1
time_steps_per_day = 24
sites = 15
time_steps = time_steps_per_day * days
months = ['01', '02', '03', '04', '05', '06', '07', '08', '09', '10', '11', '12']
years = ['2013', '2014', '2015', '2016']
#['mast number and filename prefix' ,'number of months in set' , 'number of years in set']
wind_sites_meta_data = [['wm01', 30, 4, 1 ], ['wm08', 32, 2, 1 ]]
SITES = ['SITE_1','SITE_2', 'SITE_3', 'SITE_4', 'SITE_5', 'SITE_6',
'SITE_7', 'SITE_8', 'SITE_9', 'SITE_10', 'SITE_11', 'SITE_12',
'SITE_13', 'SITE_14', 'SITE_15', 'SITE_16', 'SITE_17', 'SITE_18']
In [4]:
def demand_import(filename_or_url, import_resolution='15min',
output_resolution='15min',
day_to_year = True,
days = 365, start_date = '1/1/2015' ):#ImportDemand Data from LOCAL FILE in csv format
demand_data = pd.read_csv(filename_or_url, header = 0)
demand_data.index = pd.date_range(start_date,
periods = demand_data.iloc[:,0].size, freq = import_resolution)
if day_to_year:
demand_data = pd.concat([demand_data]*days)
demand_data.index = pd.date_range(start_date,
periods = demand_data.iloc[:,0].size, freq = import_resolution)
demand_data = demand_data.resample(output_resolution).mean()
return demand_data
# Plotly convienience functions
def make_anno1(text, fontsize, x, y):
return dict(
text=text, # annotation text
xref='paper', # use paper coordinates
yref='paper', # for both x and y coords
x=x, # x and y position
y=y, # in norm. coord.
font=go.Font(size=fontsize), # text font size
showarrow=False, # no arrow (default is True)
bgcolor='#F5F3F2', # light grey background color
bordercolor='#FFFFFF', # white borders
borderwidth=1, # border width
borderpad=fontsize # set border/text space to 1 fontsize
)
def update_axis(title):
return dict(
title=title, # axis title
tickfont=dict(size=14), # font size, default is 12
gridcolor='#FFFFFF', # white grid lines
zeroline=False, # remove thick zero line
exponentformat = 'SI'
)
def sub_dict(somedict, somekeys, default=None):
return dict([ (k, somedict.get(k, default)) for k in somekeys ])
print("---Data import and basic graphing helper functions complete---")
Data is read from .csv file into Pandas DataFrame
Cufflinks binds Pandas DataFrames directly to the plotly interactive plotting functions
In [5]:
solar_data = pd.read_csv("Data\\Solar\\SAURAN\\NMU_SOLAR.csv",
index_col = 0,
parse_dates = ['date_time'])
solar_data.head(12)
#First 12 items of the data timeseries printed from the Pandas DataFrame.
#First X rows can be seen using the head(X) function - 5 rows default
Out[5]:
In [6]:
#Select only March-2015 to March-2015 using datetime index slicing
solar_data = solar_data['03-2014':'03-2015']
#Plot timeseries using Cufflinks and Plotly directly from the Pandas DataFrame.
#Embeds fully interactive plot in the browser by simply calling the .iplot() function
solar_data["GHI"].iplot(title = 'Raw Global Horizontal and Direct Normal Irradiation from file',
yTitle = 'Irridation (W/m2)')
#Specific columns can be selected by indexing them ['ColumnName']
print("---Raw solar GHI data read succesfully and full timeseries plotted using Plotly and Cufflinks ---")
In [7]:
#Primary solar irradiation transposition and temperature loss calculations.
#Final result is kWh/kW (Hourly Capacity Factor) data
import pysolar
def filter_irradiation (irr): # Filter out all negative values and values less than 5 W/m2
if irr <= 5:
return 0
else:
return irr
lat = -34.00859
long = 25.66526
elevation = 35
tilt_angle = np.radians(34)
#Solar Zenith and Azimuth Angles calculated using Pysolar package.
#Angle taken using the midpoint of the hour (hence 30 min shift)
zenith_angles = [90 - pysolar.solar.get_altitude(lat , long,
(date+ datetime.timedelta(minutes=-30)),
elevation) for date in solar_data.index[range(solar_data.index.size)]]
azimuth_angles = [pysolar.solar.get_azimuth(lat , long,
(date+ datetime.timedelta(minutes=-30)),
elevation) for date in solar_data.index[range(solar_data.index.size)]]
solar_data['ZENITH'] = pd.Series(np.radians(zenith_angles), index = solar_data.index)
solar_data['AZIMUTH'] = pd.Series(np.radians(azimuth_angles), index = solar_data.index)
#Constants for temperature and BoP losses (soiling, wiring etc..)
other_losses = 0.95
temperature_power_coefficient = -0.0041
#Constants for transposition factors
Rr = (1 - np.sin(tilt_angle))/2
Rd = (1 + np.cos(tilt_angle))/2
rho = 0.2
#Calculate angles of incidence from hourly Zenith angles and Angles of incidence on tilted plane
solar_data['INCIDENCE'] = np.arccos( np.cos(solar_data['ZENITH'])*np.cos(tilt_angle) + \
np.sin(solar_data['ZENITH'])*np.sin(tilt_angle)*np.cos(solar_data['AZIMUTH'] - \
np.radians(180)) )
#Primary Solar irradiance combined irradiance components (Direct, Diffuse, and ground reflected)
solar_data['INCLINED'] = solar_data['DNI']*np.cos(solar_data['INCIDENCE']) + \
solar_data['DHI']*Rd + rho*solar_data['GHI']*Rr
# Filter out all negative values and values less than 5 W/m2
solar_data["INCLINED"] = solar_data["INCLINED"].apply(filter_irradiation)
#Cell temperature approximation
solar_data['TEMP_CELL'] = solar_data['TEMP'] + ((44-20)/800) * solar_data['INCLINED']
#Temperature power loss factor
solar_data['TEMP_FACTOR'] = 1 + temperature_power_coefficient*(solar_data['TEMP_CELL']-25)
#Final kWh/kW calculation hourly capacity factors for inclined panels with temperature losses
solar_data['KWH/KW'] = (solar_data['INCLINED']/1000)*solar_data['TEMP_FACTOR']*other_losses
print("--- Solar irradiance to kWh/kW computation complete---")
In [8]:
solar_data.head(12)
Out[8]:
In [9]:
solar_data[['GHI','INCLINED']].iplot(
title = 'Raw GHI, DNI, and GHI and Calculated Inclined Irradiation - Interactive Timeseries',
yTitle = 'Irradiation (W/m2)')
In [10]:
solar_data["2015"].describe()
Out[10]:
In [11]:
solar_kwh_kw_dataframe_dict = {}
for i in range(sites):
solar_kwh_kw_dataframe_dict[SITES[i]] = pd.read_csv('Data\\Solar\\MERRA-NINJA\\' +
SITES[i] +'_SOLAR.csv',
index_col = 0,
parse_dates = ['time'], # Additional index handling feilds to allow efficient date index
infer_datetime_format = True).shift(2).fillna(0) # Translate from UCT to ZA Time (+2)
print("---Processed Solar Data Files Read from csv---")
In [12]:
wind_turbine_curve = pd.read_csv('Data\Wind\TURBINE_CURVE_e400nb.csv', index_col = 0)
def select_from_curve(wind_speed):
index = wind_speed * 2
base_index = int(np.floor(index))
if base_index > 39:
base_index = 39
base_val = wind_turbine_curve.iloc[base_index][2]
kwh_kw = base_val + (index - base_index) * 0.5 * \
(wind_turbine_curve.iloc[base_index+1][2] - base_val) #interpolate between windspeed bins
return kwh_kw
print("---Wind turbine curve (Kestrel e400b) read from csv and turbine curve index function defined---")
In [13]:
print("--WASA Wind Mast Data---\n")
wind_speed_dataframe_dict = {}
sites_wasa = 1
for s in range (sites_wasa):
wind_mast = wind_sites_meta_data[s][0]
files = wind_sites_meta_data[s][1]
month = wind_sites_meta_data[s][2]
year = wind_sites_meta_data[s][3]
year_month = str("_"+years[year]+months[month])
wind_data = pd.read_csv('Data\Wind\WASA\\' + wind_mast+year_month+'_v01.csv', skiprows = 16 , index_col = 0,
parse_dates = ['date_time'],
skipinitialspace = True)[["WS_62_mean", "WS_40_mean", "WS_20_mean"]]
for i in range(files-1):
yearflag = False
month+=1
if month > 11:
month = 0
yearflag = True
if yearflag:
year+=1
yearflag = False
year_month = str("_"+years[year]+months[month])
wind_data = wind_data.append(pd.read_csv('Data\Wind\WASA\\' + wind_mast+year_month+"_v01.csv",
skiprows = 16 ,
index_col = 0,
parse_dates = ['date_time'], infer_datetime_format = True,
skipinitialspace = True)[["WS_62_mean", "WS_40_mean", "WS_20_mean"]])
#Values for logarithmic interpolation of windspeeds at different heights
z0 = 0.005
H = 20.0
H0 = 62.0
H1 = 40.0
# pad missing 62m values from previous values (minor difference if only few missing values)
# 40m and 20m values were often not recorded on many of the wind masts for varying periods of time.
# These values were logarithmically interpolated from 62m
nulls_62m = wind_data.isnull().sum()['WS_62_mean']
nulls_40m = wind_data.isnull().sum()['WS_40_mean']
nulls_20m = wind_data.isnull().sum()['WS_20_mean']
wind_data['WS_62_mean'] = wind_data['WS_62_mean'].fillna(method='pad')
wind_data['WS_40_mean'] = wind_data['WS_40_mean'].fillna(wind_data['WS_62_mean']*np.log(H1/z0)/np.log(H0/z0))
wind_data['WS_20_mean'] = wind_data['WS_20_mean'].fillna(wind_data['WS_62_mean']*np.log(H/z0)/np.log(H0/z0))
print(wind_mast+" no. of NULL values at 62m: "+str(nulls_62m))
print(wind_mast+" no. of NULL values at 40m: "+str(nulls_40m))
print(wind_mast+" no. of NULL values at 20m: "+str(nulls_20m))
print('')
wind_speed_dataframe_dict[wind_mast] = wind_data.resample('10T',
loffset = datetime.timedelta(minutes = 0)).mean()
wind_data.to_csv("Data\Wind\WASA\\" + wind_mast.upper()+'_WIND_SPEED_DATA_PROCESSED.csv')
print("---Wind speed missing data filling at 40m and 20m, logarithmically interpolated from 62m data---")
In [14]:
wind_speed_dataframe_dict['wm01']["03-2015"].iplot(
title='15min Windspeed data-set at 62, 40 and 20 meters --WASA Wind Masts --',
yTitle = 'Windspeed (m/s)')
#Only March is plotted here for more responsive interactive plotting and smaller overall notebook size
In [15]:
wind_kwh_kw_dataframe_dict = {}
sites_wasa = 1
for s in range (sites_wasa):
wind_mast = wind_sites_meta_data[s][0]
wind_data = wind_speed_dataframe_dict[wind_mast]
wind_kwh_kw_data_processed = pd.DataFrame()
wind_kwh_kw_data_processed['WS_62_kwh_kw'] = wind_data["WS_62_mean"].apply(select_from_curve).resample("H").mean()
wind_kwh_kw_data_processed['WS_40_kwh_kw'] = wind_data["WS_40_mean"].apply(select_from_curve).resample("H").mean()
wind_kwh_kw_data_processed['WS_20_kwh_kw'] = wind_data["WS_20_mean"].apply(select_from_curve).resample("H").mean()
wind_kwh_kw_dataframe_dict[wind_mast] = wind_kwh_kw_data_processed
wind_kwh_kw_data_processed.to_csv('Data\Wind\WASA\\' + wind_mast.upper()+'_WIND_KWH_KWH_DATA_PROCESSED.csv')
print("---Wind turbine power curve applied to windspeed data and converted to kWh/kW capacity factors---")
print("Sites processed successfully: ", sites_wasa)
In [21]:
wind_kwh_kw_dataframe_dict = {}
for s in range(sites_wasa):
wind_mast = wind_sites_meta_data[s][0]
wind_kwh_kw_dataframe_dict[wind_mast] = pd.read_csv('Data\Wind\WASA\\' +
wind_mast.upper()+'_WIND_KWH_KWH_DATA_PROCESSED.csv',
index_col = 0, parse_dates = ['date_time'],
infer_datetime_format = True)
print("--- Successfully imported ",sites_wasa,
" WASA wind kWh/kW hourly capacity factors from PROCESSED .csv's ---")
In [22]:
wind_kwh_kw_dataframe_dict = {}
no_to_process = 1
for s in range (no_to_process):
wind_mast = SITES[s]
wind_data = Dataset('Data\Wind\CSIR\\' + wind_mast+'.nc')
wind_array = wind_data['WSPD'][:]
np.array([row[0][0] for row in wind_array]).mean()
ind = pd.date_range(date_start, periods = 175296, freq = "15T")
wind_DataFrame = pd.DataFrame([row[3][0] for row in wind_array], index=ind, columns = ['WSPD'])
wind_data = wind_DataFrame
wind_kwh_kw_data_processed = pd.DataFrame()
#Wind data only processed to hourly capacity factors for windspeed at 40m
wind_kwh_kw_data_processed['WS_40_kwh_kw'] = wind_data["WSPD"].apply(select_from_curve).resample("H").mean()
wind_kwh_kw_dataframe_dict[wind_mast] = wind_kwh_kw_data_processed
wind_kwh_kw_data_processed.to_csv('Data\Wind\CSIR\\' +wind_mast.upper()+'_WIND_KWH_KWH_DATA_PROCESSED.csv')
print("--- Successfully Processed", no_to_process,
" site of Fraunhofer/CSIR wind speed data to kWh/kW processed hourly capacity factor .csv's ---")
In [16]:
wind_kwh_kw_dataframe_dict = {}
#Number of sites data to import
sites = 15
for i in range(sites):
wind_kwh_kw_dataframe_dict[SITES[i]] = pd.read_csv('Data\Wind\CSIR\\' +
SITES[i]+'_WIND_KWH_KWH_DATA_PROCESSED.csv',
index_col = 0,
parse_dates = True,
infer_datetime_format = True)
print("--- Successfully imported ",sites,
" Fraunhofer/CSIR wind kWh/kW processed hourly capacity factor .csv's ---")
In [17]:
# Averaged over selected timeframe ("M" - Monthly, "W" - Weekly, "D" - Daily)
interval_to_average_over = "W"
site = 0
index = wind_kwh_kw_dataframe_dict[
SITES[site]][date_start:date_end].resample(interval_to_average_over).mean().index
data = [
go.Bar(x=index,
y=solar_kwh_kw_dataframe_dict[SITES[site]]['output']
[date_start:date_end].resample(interval_to_average_over).mean(),
name='SAURAN / MERRA-2',
marker=go.Marker(color='rgb(245,160,0)')
),
go.Bar(x=index,
y=wind_kwh_kw_dataframe_dict[SITES[site]]['WS_40_kwh_kw']
[date_start:date_end].resample(interval_to_average_over).mean(),
name='WASA / CSIR',
marker=go.Marker(color='rgba(0, 70, 210, 0.96)')
)]
layout = go.Layout(
barmode='',
bargap=0.12,
#bargap=0,
plot_bgcolor='rgb(245,242,240)', # set plot and
paper_bgcolor='rgb(235,232,230)', # frame background color
height=600,
width=990,
font=go.Font(
family="Arial",
color='#635F5D',
size=13))
layout['xaxis'].update(update_axis('<b>Annual Timeseries Summed Per Selected Interval</b>'))
layout['yaxis'].update(update_axis('<b>Average Capacity Factor for Time Period</b>'))
annotations = [make_anno1('<b>Wind and Solar Capacity Factors (Averaged Weekly)</b>',
16, 0.5, 1.20)]
layout['annotations'] = annotations
fig = go.Figure(data=data, layout=layout)
po.iplot(fig)
print("---Timespan Averaged Energy Mix Plotting Complete (Monthly, Weekly, Daily, etc...) ---")
In [18]:
days = 365
demand_data_full = demand_import("Data\Scenarios\DEMAND_DATA_1_DAY.csv",
day_to_year = True, days = days, import_resolution = "H",
output_resolution='H', start_date = date_start)
demand_array = np.append(np.array(demand_data_full['Total Demand'], dtype ='float'), 0)
# Plot for single day
demand_data_full[date_start].iplot(title = 'Formulated Daily Demand Profile',
yTitle = 'Power Demand in Watts (W)')
In [19]:
litre_per_gallon = 3.78541
kwh_per_litre = 10
diesel_fuel_chart = pd.read_csv("Data\Scenarios\\DIESEL_FUEL_CHART_PROCESSED.csv")
diesel_efficiency_matrix = diesel_fuel_chart.iloc[:,[0,4,8,12,16]].as_matrix()
In [20]:
%%cython
from libc.math cimport round
cimport numpy as np
import numpy as np
import cython
from cython.parallel import prange
#Cython compiler options for speed
@cython.boundscheck(False)
@cython.nonecheck(False)
@cython.cdivision(True)
cdef hydrogen_discharge(double residual, int t,
double power_max,
double [:] hydrogen_capacity, double [:] hydrogen_discharge_power,
double hydrogen_discharge_efficiency = 0.6):
#Total hydrogen demand of power output over hour timestep
cdef double hydrogen_drain = (residual) / hydrogen_discharge_efficiency
#Set to provide power demanded and updated after logic sequence
hydrogen_discharge_power[t] = residual
#Power limitation on Fuel Cell output power
if residual > power_max:
hydrogen_drain = power_max / hydrogen_discharge_efficiency
hydrogen_discharge_power[t] = power_max
#Calculate hydrogen storage level based on energy demand compared to energy in previous timestep
cdef double residual_remainder = 0.0
hydrogen_capacity[t] = hydrogen_capacity[t-1] - hydrogen_drain
#If hydrogen was depleted or could not serve the total load power output and residual load are adjusted
if hydrogen_capacity[t] < 0:
hydrogen_discharge_power[t] = hydrogen_capacity[t-1]*hydrogen_discharge_efficiency
residual_remainder = residual - hydrogen_discharge_power[t]
hydrogen_capacity[t] = 0
return residual_remainder
@cython.boundscheck(False)
@cython.nonecheck(False)
@cython.cdivision(True)
cdef hydrogen_charge(double residual, int t,
double power_max, double max_storage_wh,
double [:] hydrogen_capacity, double [:] hydrogen_charge_power,
double hydrogen_charge_efficiency):
#power available to be used by electrolyzer in this hour
cdef double hydrogen_charging_available = -residual
#Electrolyzer power limitation
if hydrogen_charging_available > power_max:
hydrogen_charging_available = power_max
#STORAGE apacity available left for FILLING
cdef double capacity_available = max_storage_wh - hydrogen_capacity[t-1]
#Update hydrogen storage level
hydrogen_capacity[t] = hydrogen_capacity[t-1] + hydrogen_charging_available*hydrogen_charge_efficiency
#adjust averaged power output over hour if storage was filled
hydrogen_charge_power[t] = min(hydrogen_charging_available, capacity_available)
cdef double residual_remainder = residual + hydrogen_charge_power[t]
#If storage was 'over filled' or started full set to max and adjust as above
if hydrogen_capacity[t] >= max_storage_wh:
hydrogen_capacity[t] = max_storage_wh
#adjust
hydrogen_charge_power[t] = min(hydrogen_charging_available, capacity_available)
residual_remainder = residual + hydrogen_charge_power[t]
#Pass any unused energy to next level
return residual_remainder
@cython.boundscheck(False)
@cython.nonecheck(False)
@cython.cdivision(True)
cdef battery_discharge(double residual, int t,
double [:] battery_capacity, double [:] battery_discharge_power,
double battery_discharge_efficiency = 0.9592):
cdef double battery_drain = (residual) / battery_discharge_efficiency
cdef double residual_remainder = 0.0
battery_capacity[t] = battery_capacity[t-1] - battery_drain
battery_discharge_power[t] = residual
if battery_capacity[t] < 0 :
battery_capacity[t] = 0.0
battery_discharge_power[t] = battery_capacity[t-1]*battery_discharge_efficiency
residual_remainder = residual - battery_discharge_power[t]
return residual_remainder
@cython.boundscheck(False)
@cython.nonecheck(False)
@cython.cdivision(True)
cdef battery_charge(double residual, int t,
double max_storage_wh,
double [:] battery_capacity, double [:] battery_charge_power,
double battery_charge_efficiency = 0.9592):
cdef double battery_charging_available = -residual
battery_capacity[t] = battery_capacity[t-1] + \
battery_charging_available*battery_charge_efficiency
cdef double capacity_available = max_storage_wh - battery_capacity[t-1]
battery_charge_power[t] = min(battery_charging_available, capacity_available)
cdef double residual_remainder = -(battery_charging_available - battery_charge_power[t])
if battery_capacity[t] >= max_storage_wh:
battery_capacity[t] = max_storage_wh
battery_charge_power[t] = capacity_available
residual_remainder = -(battery_charging_available - (max_storage_wh - battery_capacity[t-1]))
return residual_remainder
@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
@cython.cdivision(True)
cdef diesel_run(double residual, int t,
double [:] diesel_power,
double [:] diesel_fuel_used,
double [:] diesel_efficiency_curve,
double diesel_size_kW = 4,
double minimum_load = 0.0020,
double min_size = 20,
double max_size = 2250):
cdef double residual_remainder = 0
cdef double part_load = 0
cdef double efficiency = 0
cdef int i = 0
cdef int size_index = 0
cdef int load_index = 0
if diesel_size_kW == 0:
residual_remainder = residual
else:
part_load = residual/(diesel_size_kW*1000)
diesel_power[t] = residual
# Round load factor to nearest 25% and select efficiency from matrix
load_index = <int>round(part_load*4)
efficiency = diesel_efficiency_curve[load_index]
diesel_fuel_used[t] = (residual / efficiency)/1000
if residual > diesel_size_kW*1000:
diesel_power[t] = diesel_size_kW*1000
diesel_fuel_used[t] = (diesel_power[t] / efficiency)/1000
residual_remainder = residual - diesel_power[t]
#Minimum load requirements, residual remainder to charge batteries
if part_load < minimum_load:
part_load = minimum_load # Force diesel to run at minimum load
# Round load factor to nearest 25% and select efficiency from matrix
load_index = <int>round(part_load*4)
diesel_power[t] = (part_load * diesel_size_kW)*1000
diesel_fuel_used[t] = (diesel_power[t] / efficiency)/1000
# Negative Residual Remainder for Charging from excess diesel power
residual_remainder = residual - diesel_power[t]
return residual_remainder
@cython.boundscheck(False)
@cython.nonecheck(False)
@cython.cdivision(True)
def grid_simulation(double [:] var,
double [:] demand,
double [:] solar_radiation,
double [:] wind_output,
double [:,:] diesel_efficiency_matrix,
long time_steps_per_day,
long days,
long simulation_years,
double ENS,
double hydrogen_cost,
dict tech_data,
): # Default 1hour res, No leap Year, hydrogen available
#Key Component Size Optimization Variables
cdef double pv_wp = var[0]
cdef double bat_wh = var[1]
cdef double wind_wp = var[2]
cdef double fuel_cell_wp = var[3]
cdef double electrolyzer_wp = var[4]
cdef double hydrogen_wh = var[5]
# DemandTotaling
cdef double total_period_demand = np.asarray(demand).sum()
#Diesel Sizing
cdef double diesel_size_kW = 1.2*np.asarray(demand).max()/1000
#Inverter Sizing
cdef double inverter_size_kW = 1.2*np.asarray(demand).max()/1000
#Timesteps for run
cdef int time_steps = time_steps_per_day * days
# LCOE Costing Information
# "Standard Values"
cdef double discount_rate = 0.082
cdef int lifetime = 20
cdef int n = 20
cdef double pv_cap_Wp = 2.900
cdef double pv_fom_Wp = 0.02
cdef double wind_cap_Wp = 2.720
cdef double wind_fom_Wp = 0.05
cdef double battery_cap_Wh = 0.500
cdef double battery_cap_replace_10y_Wh = 0.2500
cdef double diesel_cap_kW = 600
cdef double diesel_fuel_USD_l = 2.5
cdef double diesel_fuel_USD_kWh = 0.20
cdef double diesel_fom_kW = 15
cdef double diesel_vom_kWh = 0.015
cdef double hydrogen_storage_cap_wh = 0.055
#Seperate Hydrogen Costing variables for cost target determination
cdef double electrolyzer_cap_wp = 3.5
cdef double fuel_cell_cap_wp = 3.5
#Values from 'TECH_DATA' scanarios
discount_rate = tech_data['discount_rate']
lifetime = tech_data['lifetime']
pv_cap_Wp = tech_data['pv_cap_Wp']
pv_fom_Wp = tech_data['pv_fom_Wp']
inverter_cap_kW = tech_data['inverter_cap_kW']
inverter_fom_kW = tech_data['inverter_fom_kW']
wind_cap_Wp = tech_data['wind_cap_Wp']
wind_fom_Wp = tech_data['wind_fom_Wp']
battery_cap_Wh = tech_data['battery_cap_Wh']
battery_cap_replace_10y_Wh = tech_data['battery_cap_replace_10y_Wh']
diesel_cap_kW = tech_data['diesel_cap_kW']
diesel_fuel_USD_l = tech_data['diesel_fuel_USD_l']
diesel_fuel_USD_kWh = tech_data['diesel_fuel_USD_kWh']
diesel_fom_kW = tech_data['diesel_fom_kW']
diesel_vom_kWh = tech_data['diesel_vom_kWh']
hydrogen_storage_cap_wh = tech_data['hydrogen_storage_cap_wh']
fuel_cell_cap_wp = tech_data['fuel_cell_cap_wp']
electrolyzer_cap_wp = tech_data['electrolyzer_cap_wp']
fuel_cell_cap_wp = hydrogen_cost
electrolyzer_cap_wp = hydrogen_cost
#Costing variables derived from percentages or other variables.
cdef double battery_fom_Wh = battery_cap_Wh*0.015
cdef double electrolyzer_fom_wp = electrolyzer_cap_wp*0.05
cdef double fuel_cell_fom_wp = fuel_cell_cap_wp*0.05
cdef double minimum_diesel_load = 0.25
cdef double minimum_diesel_power = diesel_size_kW*minimum_diesel_load*1000
cdef diesel_chart_min_size = 20
cdef diesel_chart_max_size = 2250
# Component Performance Characteristics
cdef double pv_bos_eff = 0.95 # Soiling, Mismatch, Ohmic, Inter-row Shading
cdef double inverter_eff = 0.976 # Typical Inverter
cdef double battery_charge_efficiency = 0.9592 # Charging and Discharging efficiency
cdef double battery_discharge_efficiency = 0.9592
cdef double hydrogen_charge_efficiency = 0.65 # Charging and Discharging efficiency
cdef double hydrogen_discharge_efficiency = 0.6
cdef double wind_bos_eff = 1
pv_bos_eff = tech_data['pv_bos_loss_eff']
inverter_eff = tech_data['inverter_eff']
battery_charge_efficiency = tech_data['battery_charge_efficiency']
battery_discharge_efficiency = tech_data['battery_discharge_efficiency']
hydrogen_charge_efficiency = tech_data['hydrogen_charge_efficiency']
hydrogen_discharge_efficiency = tech_data['hydrogen_discharge_efficiency']
wind_bos_loss = tech_data['wind_bos_eff']
# Input Data Scaling to input Component Size
cdef double solar_array_size_wp = pv_wp
cdef double [:] solar_generation = np.asarray(solar_radiation) * solar_array_size_wp * pv_bos_eff
cdef double solar_energy_sum = np.asarray(solar_generation).sum()
cdef double wind_array_size_wp = wind_wp
cdef double [:] wind_generation = np.asarray(wind_output) * wind_array_size_wp * wind_bos_eff
cdef double wind_energy_sum = np.asarray(wind_generation).sum()
# Initialize Energy Time Vectors
cdef double [:] battery_capacity = np.zeros(time_steps + 1)
cdef double [:] battery_discharge_power = np.zeros(time_steps + 1)
cdef double [:] battery_charge_power = np.zeros(time_steps + 1)
cdef double [:] hydrogen_capacity = np.zeros(time_steps + 1)
cdef double [:] hydrogen_discharge_power = np.zeros(time_steps + 1)
cdef double [:] hydrogen_charge_power = np.zeros(time_steps + 1)
cdef double [:] electrolyser_power = np.zeros(time_steps + 1)
cdef double [:] fuel_cell_power = np.zeros(time_steps + 1)
cdef double [:] diesel_power = np.zeros(time_steps + 1)
cdef double [:] diesel_fuel_used = np.zeros(time_steps + 1)
cdef double [:] solar_generation_utilized = np.zeros(time_steps + 1)
cdef double [:] wind_generation_utilized = np.zeros(time_steps + 1)
cdef double [:] renewable_generation_curtailed = np.zeros(time_steps + 1)
cdef double [:] energy_not_served = np.zeros(time_steps + 1)
cdef double [:] residual_load = np.zeros(time_steps + 1)
# battery_capacity Max Cap and starting value
cdef double max_storage_wh = bat_wh
cdef double max_hydrogen_storage_wh = hydrogen_wh
battery_capacity[-1] = max_storage_wh/2
hydrogen_capacity[-1] = max_hydrogen_storage_wh/2
# annulisation of capital expenditure
def annualise(double pv, double i, int n):
return <double>(pv * ((i * (1 + i) ** n) / ((1 + i) ** n - 1)))
cdef int i = 0
# EFFICIENCY CURVE CALCULATION & FUEL USE
# if generator size is less than or greater than max set to minimum or maximum
if diesel_size_kW < diesel_chart_min_size:
diesel_size_eff = diesel_chart_min_size
if diesel_size_kW > diesel_chart_max_size:
diesel_size_eff = diesel_chart_max_size
for i in range(24): # Find disel generator size from matrix
if diesel_size_eff >= diesel_efficiency_matrix[i][0]:
size_index = i
break
diesel_efficiency_curve = diesel_efficiency_matrix[size_index]
cdef double [:] renewable_generation = (np.asarray(solar_generation) + \
np.asarray(wind_generation))*inverter_eff
residual_load = np.asarray(demand) - np.asarray(renewable_generation)
cdef int t = 0
cdef double Residual, Demand, residual_remainder = 0
cdef double battery_drain = 0
cdef double battery_charging_available = 0
cdef double capacity_available = 0
cdef double diesel_demand = 0
cdef double fuel_cell_power_max = fuel_cell_wp
cdef double electrolyzer_power_max = electrolyzer_wp
cdef hydrogen_energy_available = 0
# Calculate Mini-Grid Operation for each timestep
for t in range(time_steps):
# Temp Residual for speed
Residual = residual_load[t]
Demand = demand[t]
residual_remainder = 0.0
hydrogen_energy_available = min(hydrogen_capacity[t-1]*hydrogen_discharge_efficiency,
fuel_cell_power_max)
capacity_available = battery_capacity[t-1]*battery_discharge_efficiency + \
hydrogen_energy_available
if Residual > 0: # Discharging through battery and inverter or served by diesel
# Load Served directly by solar PV & wind
solar_generation_utilized[t] = solar_generation[t] * inverter_eff
wind_generation_utilized[t] = wind_generation[t] * inverter_eff
if Residual > capacity_available: #Diesel will need to run
if Residual < minimum_diesel_power:
residual_remainder = diesel_run(Residual, t,
diesel_power, diesel_fuel_used,
diesel_efficiency_curve,
diesel_size_kW,
minimum_diesel_load)
residual_remainder = battery_charge(residual_remainder, t, max_storage_wh,
battery_capacity,
battery_charge_power,
battery_charge_efficiency)
residual_remainder = hydrogen_charge(residual_remainder, t,
electrolyzer_power_max, max_hydrogen_storage_wh,
hydrogen_capacity, hydrogen_charge_power,
hydrogen_charge_efficiency)
elif Residual > minimum_diesel_power:
residual_remainder = battery_discharge((Residual - minimum_diesel_power), t,
battery_capacity, battery_discharge_power)
residual_remainder = hydrogen_discharge(residual_remainder, t, fuel_cell_power_max,
hydrogen_capacity, hydrogen_discharge_power,
hydrogen_discharge_efficiency)
residual_remainder = diesel_run((Residual - battery_discharge_power[t] - \
hydrogen_discharge_power[t]), t,
diesel_power, diesel_fuel_used,
diesel_efficiency_curve,
diesel_size_kW,
minimum_diesel_load)
energy_not_served[t] = residual_remainder
else: # Storage can meet demand
residual_remainder = battery_discharge(Residual, t,
battery_capacity,
battery_discharge_power,
battery_discharge_efficiency)
residual_remainder = hydrogen_discharge(residual_remainder, t, fuel_cell_power_max,
hydrogen_capacity, hydrogen_discharge_power,
hydrogen_discharge_efficiency)
# Negative Residual: Serve load and charge battery with remainder bypassing inverter loss
if Residual < 0:
if solar_generation[t] > Demand:
solar_generation_utilized[t] = Demand
elif solar_generation[t] < Demand:
solar_generation_utilized[t] = solar_generation[t]*inverter_eff
wind_generation_utilized[t] = Demand - solar_generation_utilized[t]
residual_remainder = battery_charge(Residual, t, max_storage_wh,
battery_capacity,
battery_charge_power,
battery_charge_efficiency)
residual_remainder = hydrogen_charge(residual_remainder, t,
electrolyzer_power_max, max_hydrogen_storage_wh,
hydrogen_capacity, hydrogen_charge_power,
hydrogen_charge_efficiency)
renewable_generation_curtailed[t] = -residual_remainder
if residual_remainder < -0.1:
if (solar_generation[t] > Demand) & (wind_generation[t] > Demand):
solar_generation_utilized[t] = Demand/2
wind_generation_utilized[t] = Demand/2
if (solar_generation[t] > Demand) & (wind_generation[t] < Demand):
solar_generation_utilized[t] = Demand - wind_generation[t]/2
wind_generation_utilized[t] = wind_generation[t]/2
if (solar_generation[t] < Demand) & (wind_generation[t] > Demand):
solar_generation_utilized[t] = solar_generation[t]/2
wind_generation_utilized[t] = Demand - solar_generation[t]/2
if (solar_generation[t] < Demand) & (wind_generation[t] < Demand):
solar_generation_utilized[t] = Demand/2
wind_generation_utilized[t] = Demand/2
# Sum Relevant Timeseries as Indicators
cdef double renewable_curtailed_total = np.asarray(renewable_generation_curtailed).sum()
cdef double diesel_energy_total = np.asarray(diesel_power).sum()
cdef double battery_energy_total = np.asarray(battery_discharge_power).sum()
cdef double hydrogen_energy_total = np.asarray(hydrogen_discharge_power).sum()
cdef double solar_direct_energy_total = np.asarray(solar_generation_utilized).sum()
if solar_direct_energy_total == 0:
solar_direct_energy_total = 1
cdef double wind_direct_energy_total = np.asarray(wind_generation_utilized).sum()
if wind_direct_energy_total == 0:
wind_direct_energy_total = 1
cdef double renewable_energy_total = battery_energy_total + hydrogen_energy_total + \
solar_direct_energy_total + wind_direct_energy_total + \
hydrogen_energy_total
cdef double revenue_energy_total = renewable_energy_total + diesel_energy_total
cdef double diesel_fuel_total = np.asarray(diesel_fuel_used).sum()
if diesel_fuel_total == 0:
diesel_fuel_total = 1
if diesel_fuel_total < 0:
diesel_fuel_total = 1
cdef double prod_total = revenue_energy_total / 1000 / 100
cdef double diesel_energy_percent = diesel_energy_total / total_period_demand * 100
cdef double total_curtailed_percent = (renewable_curtailed_total / renewable_energy_total) * 100
cdef double diesel_average_efficiency = (diesel_energy_total / 1000) / diesel_fuel_total
cdef double lc_ENS = 0
if diesel_energy_percent > ENS:
lc_ENS = 100000
# Individual Levelised Components in $/kWh
cdef double lc_pv_cap = annualise(solar_array_size_wp * pv_cap_Wp,
discount_rate, lifetime) / prod_total * simulation_years
cdef double lc_pv_om = pv_fom_Wp * solar_array_size_wp / prod_total *simulation_years
cdef double lc_inverter_cap = annualise(inverter_cap_kW * inverter_size_kW,
discount_rate, lifetime) / prod_total * simulation_years
cdef double lc_inverter_fom = inverter_size_kW * inverter_fom_kW / prod_total * simulation_years
cdef double lc_wind_cap = annualise(wind_array_size_wp * wind_cap_Wp,
discount_rate, lifetime) / prod_total * simulation_years
cdef double lc_wind_om = wind_fom_Wp * wind_array_size_wp / prod_total * simulation_years
cdef double lc_battery_cap = annualise(max_storage_wh * (battery_cap_Wh + battery_cap_replace_10y_Wh),
discount_rate, lifetime) / prod_total * simulation_years
cdef double lc_battery_om = battery_fom_Wh * max_storage_wh / prod_total * simulation_years
cdef double lc_diesel_cap = annualise(diesel_cap_kW * diesel_size_kW,
discount_rate, lifetime) / prod_total * simulation_years
cdef double lc_diesel_om = ((diesel_fom_kW * diesel_size_kW) + (diesel_vom_kWh * diesel_energy_total / 1000)) / prod_total * simulation_years
cdef double lc_diesel_fuel = diesel_fuel_USD_kWh * diesel_fuel_total / prod_total * simulation_years
cdef double lc_fuel_cell_cap = annualise(fuel_cell_wp * fuel_cell_cap_wp,
discount_rate, lifetime) / prod_total * simulation_years
cdef double lc_fuel_cell_om = fuel_cell_wp * fuel_cell_fom_wp / prod_total * simulation_years
cdef double lc_electrolyzer_cap = annualise(electrolyzer_wp * electrolyzer_cap_wp,
discount_rate, lifetime) / prod_total * simulation_years
cdef double lc_electrolyzer_om = electrolyzer_wp * electrolyzer_fom_wp / prod_total * simulation_years
cdef double lc_hydrogen_storage_cap = annualise(hydrogen_storage_cap_wh * hydrogen_wh,
discount_rate, lifetime) / prod_total * simulation_years
# LCOE Total in $c/kWh
cdef double LCOE = ((lc_pv_cap + lc_pv_om) + \
(lc_wind_cap + lc_wind_om) + \
(lc_inverter_cap + lc_inverter_fom) + \
(lc_battery_cap + lc_battery_om) + \
(lc_diesel_cap + lc_diesel_om + lc_diesel_fuel) + \
(lc_fuel_cell_cap + lc_fuel_cell_om + \
lc_electrolyzer_cap + lc_electrolyzer_om + \
lc_hydrogen_storage_cap + \
(lc_ENS)))
# Returns Dictionary of all relevant indicators and Timeseries
return {"LCOE": LCOE,
"demand": np.asarray(demand),
"residual_load": np.asarray(residual_load),
"solar_generation": np.asarray(solar_generation),
"solar_generation_utilized": np.asarray(solar_generation_utilized),
"wind_generation" : np.asarray(wind_generation),
"wind_generation_utilized" : np.asarray(wind_generation_utilized),
"renewable_generation_curtailed": np.asarray(renewable_generation_curtailed),
"renewable_generation_utilized" : renewable_energy_total,
"battery_capacity": np.asarray(battery_capacity),
"battery_charge_power": np.asarray(battery_charge_power),
"battery_discharge_power": np.asarray(battery_discharge_power),
"hydrogen_capacity": np.asarray(hydrogen_capacity),
"hydrogen_charge_power": np.asarray(hydrogen_charge_power),
"hydrogen_discharge_power": np.asarray(hydrogen_discharge_power),
"diesel_power": np.asarray(diesel_power),
"diesel_fuel_used": np.asarray(diesel_fuel_used),
"diesel_energy_percent": diesel_energy_percent,
"diesel_average_efficiency": diesel_average_efficiency,
"total_curtailed_percent": total_curtailed_percent,
"LCOE_components": {"lc_pv_cap": lc_pv_cap,
"lc_pv_om": lc_pv_om,
"lc_wind_cap" : lc_wind_cap,
"lc_wind_om" : lc_wind_om,
"lc_inverter_cap" : lc_inverter_cap,
"lc_inverter_fom" : lc_inverter_fom,
"lc_battery_cap": lc_battery_cap,
"lc_battery_om": lc_battery_om,
"lc_diesel_cap": lc_diesel_cap,
"lc_diesel_om": lc_diesel_om,
"lc_diesel_fuel": lc_diesel_fuel,
"lc_fuel_cell_cap": lc_fuel_cell_cap,
"lc_fuel_cell_om": lc_fuel_cell_om,
"lc_electrolyzer_cap" : lc_electrolyzer_cap,
"lc_electrolyzer_om" : lc_electrolyzer_om,
"lc_hydrogen_storage_cap" : lc_hydrogen_storage_cap
}
}
print("---Main mini-grid annual simulation and costing function compiled into C using Cython---")
In [23]:
#Random particle position generator function with veolocity limits and dimensioning
def init_X(n, dim, v_range):
X = np.random.uniform(1,v_range, size=(n,dim))
return X
def pso(function, dim, particles, generations, demand, solar_radiation,
wind_output, diesel_efficiency_matrix,
ranges, time_steps_per_day, days, simulation_years, ENS, hydrogen_cost,
tech_data):
#PSO Constants for social/cognitive weights, inertial weights, velocity caps.
#Typical "Standard" values for simple PSO
c1 = 1.49445
c2 = 1.49445
wstart = 1.2
wend = 0.7
w_vec = np.ones(shape=(particles, dim))
w_min = 0.6
w_max = 2.7
v_max = ranges * 0.15 # Vmax set to 15% of dynamic range
#Initial Random Particle Position Generation
X = init_X(particles, dim, ranges)
#additional random points within ranges used to set initial velocity vectors towards
Xt = init_X(particles, dim, ranges)
V = (X-Xt)/2 #Initialize Velocity towards random points and half vector
p_best = np.ones(shape=(particles, dim))*1000
l_best = np.ones(shape=(particles, dim))*1000
p_best_val = np.ones(particles)*1000
g_best = np.array([0.0, 0.0, 0.0, 0.0, 0.0, 0.0])
g_best_val = 10000
generation_bests = np.zeros(generations)
#Chaos random number seeds
r1 = 0.1123
r2 = 0.3123
for iteration in range(generations):
if iteration > 0:
#X_previous = X # Keep previous Positions for bound handling techniques
X+= V
#FIX PARTICLES THAT GO OUT OF BOUNDS
for i in range(particles):
for n in range(dim):
#if bottom of the range is hit, set to 0
if X[i][n] < 0:
X[i][n] = 0
V[i][n] = X[i][n] - g_best[n] #with velocity accelerated toward best particle
#if top of the range is exceeded, reinitialize offending position component
#with random velocity direction,
if X[i][n] > ranges:
X[i][n] = np.random.uniform(0,ranges)
V[i][n] = X[i][n] - g_best[n] #with velocity accelerated toward best particle
#temporary particle storage array
results_temp = np.zeros(particles)
for p in range(particles):
results_temp[p] = grid_simulation(X[p], demand_array, solar_radiation, wind_output,
diesel_efficiency_matrix,
time_steps_per_day, days,
simulation_years, ENS, hydrogen_cost,
tech_data)["LCOE"]
#Average and minimum fitness of particles caluclated for Dynamic Intertial Weight Adjustment
#Calculate average and minimum values to adjust indivitual particles inertial weights.
#Adjust each particle's inertial weights based on it's ranking/comparison to the global average.
#If particles are below the average fitness set inertia high to explore further,
#if above average set inertia low for closer exploititive search
ave = results_temp.mean()
min = results_temp.min()
#Main particle set evaluation loop
for i in range(particles):
#retrieve i'th particle
res = results_temp[i]
#Calculate and set dynamic inertial weights
narrow = w_min + (((w_max - w_min)*(res - min))/(ave - min))
wide = w_max
if res <= ave:
w_vec[i] = np.array([narrow]*dim)
else:
w_vec[i] = np.array([wide]*dim)
#Compare particle current fitness to particles own best and replace/store value if better
if res < p_best_val[i]:
p_best[i] = np.array(X[i])
p_best_val[i] = res
#Compare, store and replace for global best particle position and value
if res < g_best_val:
g_best = np.array(X[i])
g_best_val = res
g_best_index = i
#Chaos random variables r1 and r2
r1 = 4*r1*(1-r1)
r2 = 4*r2*(1-r2)
#Velocity Vector Update formula with dynamic intertial weight and chaos random numbers
V = (w_vec* V) + (c1 * r1 * (p_best - X)) + (c2 * r2 * (g_best - X))
#PARTICLE VELOCITY LIMITATION checking to v_max
#Total vector velocity is calculated and capped, not individual variables
Vlen = []
for k in range(particles):
vel = np.sqrt(V[k].dot(V[k]))
Vlen.append(vel)
if vel > v_max:
V[k] *= v_max / vel
#Output for final Best solution found with variable and objective function values
print('Best particle: ', g_best, 'with value of: %0.4f' % g_best_val)
#Output format of tech sizes array from PSO is:
#[Solar PV (Wp), Battery (Wh), Wind(Wp), Fuel Cell (Wp), Electrolyser (Wp) H2 Storage (Wh)]
#Return best value and variables
return (g_best, g_best_val)
print('---Particle Swarm Optimization Function Successfully Defined (Standard Python: No C code)---')
In [25]:
tech_data_DataFrame = pd.read_csv('Data\Scenarios\TECH_DATA_ONE.csv')
tech_data_dict = tech_data_DataFrame.T.to_dict()
SCENARIOS = tech_data_DataFrame.index
tech_data_DataFrame
Out[25]:
In [27]:
cost_steps = np.linspace(10, 0.25, 11)
scenarios = tech_data_DataFrame.index.size
sites = tech_data_DataFrame['site'].unique().size
cost_target_iterations = 1
site_results_List = [0]*scenarios
scenario_results_List = [0]*scenarios
ENS = 200
show_sizes_output = False
print('----------- Particle Swarm Optimization Runs Started -----------')
#Loop through all modelled sites
for s in range(scenarios):
SITE = tech_data_DataFrame.iloc[s]['site']
TECHNOLOGIES = tech_data_DataFrame.iloc[s]['technologies']
YEAR = tech_data_DataFrame.iloc[s]['year']
# read solar and wind kwh/kw data from processed files and input into numpy arrays
solar = np.append(np.asarray(solar_kwh_kw_dataframe_dict[SITE]["output"]
[date_start : date_end]), 0)
wind = np.append(np.asarray(wind_kwh_kw_dataframe_dict[SITE] ['WS_40_kwh_kw']
[date_start : date_end]), 0)
print('--- SITE: '+str(SITE)+' | TECHNOLOGIES : '+str(TECHNOLOGIES)+' | YEAR: '+str(YEAR) +' ---')
H2Perc = 0
FuelCellSize = 0
#Run optimization algorithm 3 times for each scenario and take the best result
#Also allows for evaluation of inter-run PSO performance
runs = 2
optimal_test_results = []
run_traces = []
#number of particles and generations used in runs
particles = 500
generations = 50
#Shorten diesel only simulation significantly
if TECHNOLOGIES== 'D':
particles = 20
generations = 10
ranges = 600000 #dynamic range of search region - (set annectotally here at the moment)
hydrogen = 1 # Option to include hydrogen as an option or not
hydrogen_cost = cost_steps[0] # Hydrogen cost target steps
# print('----------- Hydrogen Cost: '+ str(hydrogen_cost) + '-----------')
min_lcoe = 1000 # High initial value to compare run results against to find minimum
min_index = 0
tech_data = tech_data_dict[s]
for i in range(runs):
#Main call to Optimization Routine
x1 = pso(grid_simulation, 6, particles, generations,
demand_array, solar, wind,
diesel_efficiency_matrix, ranges,
time_steps_per_day, days,
simulation_years, ENS, hydrogen_cost,
tech_data)
#List to collate multiple PSO runs
optimal_test_results.append(x1)
#Compare and save index of best result
if optimal_test_results[i][1] < min_lcoe:
min_lcoe = optimal_test_results[i][1]
min_index = i
#Resimulation of best grid result to capture outputs into dataframes for analysis & plotting etc
run_result = grid_simulation(optimal_test_results[min_index][0],
demand_array, solar, wind,
diesel_efficiency_matrix, 24, 365, 1, ENS, hydrogen_cost,
tech_data)
#Shorten diesel only simulation significantly
if TECHNOLOGIES == 'D':
run_result = grid_simulation(np.array([0.0, 0.0, 0.0, 0.0, 0.0, 0.0]),
demand_array, solar, wind,
diesel_efficiency_matrix, 24, 365, 1, ENS, hydrogen_cost,
tech_data)
index = pd.date_range(date_start, periods = time_steps+1, freq = "H")
df = pd.DataFrame(sub_dict(run_result,
('demand', 'solar_generation',
'wind_generation', 'residual_load',
'solar_generation_utilized', 'wind_generation_utilized',
'renewable_generation_curtailed', 'battery_charge_power',
'battery_discharge_power', 'battery_capacity',
'hydrogen_charge_power', 'hydrogen_discharge_power',
'hydrogen_capacity', 'diesel_power')
),
index = index,
)
#Capture optimal grid sizing configuration in simple array
tech_size_array = optimal_test_results[min_index][0]
tech_size_array = np.append(tech_size_array , run_result['LCOE'])
# Renewable percent 100-diesel%
tech_size_array = np.append(tech_size_array ,
100 - 100*df['diesel_power'].sum()/df['demand'].sum())
#Energy used 100-curatiled%
tech_size_array = np.append(tech_size_array ,
100 - 100*df['renewable_generation_curtailed'].sum()/df['demand'].sum())
#Hydrogen %
tech_size_array = np.append(tech_size_array ,
100*df['hydrogen_discharge_power'].sum()/df['demand'].sum())
#Hydrogen Cost
tech_size_array = np.append(tech_size_array,
hydrogen_cost)
scenario_results_List[s] = tech_size_array
H2Perc = 100*df['hydrogen_discharge_power'].sum()/df['demand'].sum()
DiesPerc = 100*df['diesel_power'].sum()/df['demand'].sum()
CurtPerc = 100*df['renewable_generation_curtailed'].sum()/df['demand'].sum()
FuelCellSize = tech_size_array[3]
#Simple Output Results printing
print("%.4f" % min_lcoe, ' c$/kWh - From run number ', min_index+1)
if show_sizes_output: # Can suppress output for larger sets of scenarios
for i in range(runs):
print('run ', i+1 ,'result: ', optimal_test_results[i][1])
print('')
print('PV size : ', "%.2f" % (optimal_test_results[min_index][0][0]/1000), ' kWp')
print('Battery size : ', "%.2f" % (optimal_test_results[min_index][0][1]/1000), ' kWh')
print('Wind size : ', "%.2f" % (optimal_test_results[min_index][0][2]/1000), ' kWp')
print('FC size : ', optimal_test_results[min_index][0][3])
print('Ely size : ', optimal_test_results[min_index][0][4])
print('H2 Store size: ', optimal_test_results[min_index][0][5])
print('H2 % : ',
100*df['hydrogen_discharge_power'].sum()/df['demand'].sum())
print('Renewable % : ', "%.2f" %
(100 - 100*df['diesel_power'].sum()/df['demand'].sum()), ' %')
print('Curtailed % : ', "%.2f" %
(100*df['renewable_generation_curtailed'].sum()/df['demand'].sum()), ' %')
if (FuelCellSize >= 2500) or (H2Perc >= 10):
print("-----:::::::---- Hydrogen Economic Breakeven: ",
hydrogen_cost,"$/kW ---::::::::-----")
print("")
print("")
break
print("---Mini-Grid Sizing Computation Complete---\n")
#Output format of tech sizes array from PSO is:
#[Solar PV (Wp), Battery (Wh), Wind(Wp), Fuel Cell (Wp), Electrolyser (Wp) H2 Storage (Wh)]
In [28]:
# 7.3 Output variables list for column headings of results DataFrame output for writing to file
OUTPUT_VARIABLES = ['Site', 'Technologies', 'Year',
'PV Size', 'Battery Size', 'Wind Size', 'Fuel Cell Size', 'Electrolyser Size',
'Hydrogen Store Size', 'LCOE', 'Renewable Penetration %', 'Generated Energy Used %',
'Hydrogen Share %',"Hyrdrogen Breakeven Cost",
'PV Cost', 'Wind Cost', 'Battery Cost', 'Fuel Cell Cost',
'Electrolyzer Cost','H2 Store', 'Diesel Cost',
'PV Energy', 'Wind Energy', 'Battery Energy', 'Hydrogen Energy',
'Diesel Energy', 'Curtailed Energy',
'Total Demand']
In [29]:
OUPUT_FILENAME = "SINGLE_CSV_RESULTS.csv"
scenarios = tech_data_DataFrame.index.size
sites = tech_data_DataFrame['site'].unique().size
cost_target_iterations = 1
ENS = 200
show_sizes_output = False
scenario_results_List_2 = scenario_results_List.copy() #prevent acciental overwrite - temp measure
writing_results_records_List = [0]*scenarios
print('----------- Writing Optimization Runs to CSV file on disk -----------')
for s in range(scenarios):
SITE = tech_data_DataFrame.iloc[s]['site']
TECHNOLOGIES = tech_data_DataFrame.iloc[s]['technologies']
YEAR = tech_data_DataFrame.iloc[s]['year']
solar = np.append(np.asarray(solar_kwh_kw_dataframe_dict[SITE]["output"]
[date_start : date_end]), 0)
wind = np.append(np.asarray(wind_kwh_kw_dataframe_dict[SITE]['WS_40_kwh_kw']
[date_start : date_end]), 0)
tech_sizes = scenario_results_List_2[s]
hydrogen_cost = tech_sizes[10]
run_result = grid_simulation(tech_sizes,
demand_array, solar, wind,
diesel_efficiency_matrix, 24, 365, 1, ENS , hydrogen_cost,
tech_data)
index = pd.date_range(date_start, periods = time_steps+1, freq = "H")
df = pd.DataFrame(sub_dict(run_result,
('demand' , 'solar_generation', 'wind_generation', 'residual_load',
'solar_generation_utilized', 'wind_generation_utilized',
'renewable_generation_curtailed', 'battery_charge_power',
'battery_discharge_power', 'battery_capacity',
'hydrogen_charge_power', 'hydrogen_discharge_power',
'hydrogen_capacity' , 'diesel_power')
),
index = index)
LCOE_components = run_result['LCOE_components']
LCOE_solar_total = LCOE_components['lc_pv_cap'] + LCOE_components['lc_pv_om']
LCOE_wind_total = LCOE_components['lc_wind_cap'] + LCOE_components['lc_wind_om']
LCOE_battery_total = LCOE_components['lc_battery_cap'] + LCOE_components['lc_battery_om']
LCOE_hydrogen_fuel_cell_total = LCOE_components['lc_fuel_cell_cap'] + \
LCOE_components['lc_fuel_cell_om']
LCOE_hydrogen_electrolyzer_total = LCOE_components['lc_electrolyzer_cap'] + \
LCOE_components['lc_electrolyzer_om']
LCOE_hydrogen_storage_total = LCOE_components['lc_hydrogen_storage_cap']
LCOE_hydrogen_total = LCOE_hydrogen_fuel_cell_total + \
LCOE_hydrogen_electrolyzer_total + \
LCOE_hydrogen_storage_total
LCOE_diesel_total = LCOE_components['lc_diesel_cap'] + \
LCOE_components['lc_diesel_om'] + \
LCOE_components['lc_diesel_fuel']
scenario_results_List_temp = scenario_results_List_2[s].copy().tolist()
#inserted at start of list (in reverse order) to align columns at the start of the DataFrame/CSV Table
scenario_results_List_temp.insert(0, YEAR)
scenario_results_List_temp.insert(0, TECHNOLOGIES)
scenario_results_List_temp.insert(0, SITE)
scenario_results_List_temp.append(LCOE_solar_total)
scenario_results_List_temp.append(LCOE_wind_total)
scenario_results_List_temp.append(LCOE_battery_total)
scenario_results_List_temp.append(LCOE_hydrogen_fuel_cell_total)
scenario_results_List_temp.append(LCOE_hydrogen_electrolyzer_total)
scenario_results_List_temp.append(LCOE_hydrogen_storage_total)
scenario_results_List_temp.append(LCOE_diesel_total)
#Energy Share percentages
demand = df['demand'].sum()
scenario_results_List_temp.append(df['solar_generation_utilized'].sum())
scenario_results_List_temp.append(df['wind_generation_utilized'].sum())
scenario_results_List_temp.append(df['battery_discharge_power'].sum())
scenario_results_List_temp.append(df['hydrogen_discharge_power'].sum())
scenario_results_List_temp.append(df['diesel_power'].sum())
scenario_results_List_temp.append(df['renewable_generation_curtailed'].sum())
scenario_results_List_temp.append(demand)
writing_results_records_List[s] = scenario_results_List_temp
print('--- SITE: '+str(SITE)+' | TECHNOLOGIES : '+str(TECHNOLOGIES)+' | YEAR: '+str(YEAR) +' --- WRITTING')
#Save to CSV file
res = pd.DataFrame(writing_results_records_List, columns = OUTPUT_VARIABLES)
res.to_csv("Results\\" + OUPUT_FILENAME)
print('----------- Writing Optimization Runs to CSV file on disk ----------- COMPLETE')
print('-----------', scenarios, ' runs written to disk')
In [34]:
plot_grid_timeseries(date_start, 15, df)
In [35]:
interval = "W"
plot_energy_use_profile(tech_sizes, interval, True)
In [36]:
interval = "D"
plot_energy_use_profile(tech_sizes, interval, True)
In [37]:
plot_LCOE_components(run_result['LCOE_components'])
In [31]:
#IF PERFORMANCE IS POOR REDUCE DURATION, REMOVE UNECCESARY TIMESERIES
#OR PLOT IN STANDALONE WINDOW po.plot() instead of po.iplot()
def plot_grid_timeseries(start_date, days, data_frame):
plot_days = days
plot_time_steps = time_steps_per_day * plot_days
index = pd.date_range(date_start, periods = plot_time_steps, freq = "H")
data = [
go.Scatter(x=index,
y = data_frame['demand'],
name='Demand - Total',
line=go.Line(color='rgba(0, 0, 0, 1)',
dash='',
width=4.75,
shape='hvh')),
go.Bar(x=index,
y = data_frame['solar_generation_utilized'],
name='PV - Direct',
marker=go.Marker(color='rgb(235,150,0)')),
go.Bar(x=index,
y= data_frame['wind_generation_utilized'],
name='Wind - Direct',
marker=go.Marker(color='rgba(0, 70, 210, 0.96)')),
go.Scatter(x=index,
y= data_frame['solar_generation'],
name='PV - Total',
line=go.Line(color='gold',
width=4.5,
shape='spline')),
go.Scatter(x=index,
y= data_frame['solar_generation']+df['wind_generation'],
name='PV & Wind - Total',
line=go.Line(color='rgba(80, 0, 170, 1)',
dash='dash',
width=4.7,
shape='spline')),
go.Scatter(x=index,
y= data_frame ['wind_generation'],
name='Wind - Total',
line=go.Line(color='rgba(153, 200, 210, 1)',
width=4.5,
shape='spline')),
go.Scatter(x=index,
y=data_frame['battery_capacity'],
name='Battery - SOC',
line=go.Line(color='rgba(5, 15, 120, 1)',
dash ="dashdot",
width=5.3,
shape='spline'),
yaxis = 'y2'),
go.Scatter(x=index,
y=df['hydrogen_capacity'],
name='H2 Level',
line=go.Line(color='rgba(245, 45, 42, 1)',
dash ="dashdot",
width=5.3,
shape='spline'),
yaxis = 'y2'),
go.Bar(x=index,
y=data_frame['battery_discharge_power'],
name='Batt - Discharge',
marker=go.Marker(color='rgba(150, 200, 0, 1)')),
go.Bar(x=index,
y=df['hydrogen_discharge_power'],
name='Fuel Cell',
marker=go.Marker(color='rgba(255, 100, 100, 1)')),
go.Bar(x=index,
y=data_frame['diesel_power'],
name='Diesel - Power',
marker=go.Marker(color='rgba(100, 45, 30, 1)')),
go.Bar(x=index,
y=data_frame['battery_charge_power'],
name='Batt - Charge',
marker=go.Marker(color='rgba(0, 75, 0, 1)')),
go.Bar(x=index,
y=data_frame['hydrogen_charge_power'],
name='H2 - Charge',
marker=go.Marker(color='rgba(0, 185, 110, 1)')),
go.Bar(x=index,
y=data_frame["renewable_generation_curtailed"],
name='RE - Curtailed',
marker=go.Marker(color='rgb(230, 100, 0)')
)]
layout = go.Layout(
barmode='stack',
bargap=0.15,
plot_bgcolor='rgb(245,242,240)', # plot area color
paper_bgcolor='rgb(235,232,230)', # frame background color
height=550,
width=800,
font=go.Font(
family="Arial",
color='black',
size=14),
yaxis=dict(),
yaxis2=dict(overlaying='y', side='right',))
layout['xaxis'].update(update_axis('<b>Date-Time : Hourly Resolution</b>'))
layout['yaxis'].update(update_axis('<b>Power (W)</b>'))
annotations = [make_anno1('<b>Mini-Grid Temporal Energy Simulation Timeseries</b>',
16, 0.5, 1.23)]
layout['annotations'] = annotations
fig = go.Figure(data=data, layout=layout)
po.iplot(fig)
In [32]:
def plot_energy_use_profile(tech_sizes, interval_to_average_over, display_interactive):
run_result = grid_simulation(tech_sizes,
demand_array, solar, wind,
diesel_efficiency_matrix, 24, 365, 1, ENS, hydrogen_cost,
tech_data)
index = pd.date_range(date_start, periods = time_steps+1, freq = "H")
df = pd.DataFrame(sub_dict(run_result,
('demand', 'solar_generation', 'wind_generation', 'residual_load',
'solar_generation_utilized', 'wind_generation_utilized',
'renewable_generation_curtailed', 'battery_charge_power',
'battery_discharge_power', 'battery_capacity',
'hydrogen_charge_power', 'hydrogen_discharge_power',
'hydrogen_capacity' , 'diesel_power')
),
index = index)
index = pd.date_range(date_start, periods = time_steps+1, freq = "H")
df_annual_timeseries = pd.DataFrame(
sub_dict(
run_result, ('solar_generation_utilized', 'wind_generation_utilized',
'battery_discharge_power', 'hydrogen_discharge_power',
'diesel_power', 'renewable_generation_curtailed')),
index = index)
df_monthly_sums = df_annual_timeseries[0:-1].resample(interval_to_average_over).sum()/1000
index = df_monthly_sums.index
cf.set_config_file(theme='pearl')
data = [
go.Bar(x=index, y=df_monthly_sums['solar_generation_utilized'],
name='Solar',
marker=go.Marker(color='rgb(245,160,0)')),
go.Bar(x=index, y=df_monthly_sums['wind_generation_utilized'],
name='Wind',
marker=go.Marker(color='rgba(0, 70, 210, 0.96)')),
go.Bar(x=index, y=df_monthly_sums['battery_discharge_power'],
name='Battery',
marker=go.Marker(color='rgba(110, 185, 30, 1)')),
go.Bar(x=index, y=df_monthly_sums['hydrogen_discharge_power'],
name='Hydrogen',
marker=go.Marker(color='rgba(255, 100, 100, 1)')),
go.Bar(x=index, y=df_monthly_sums['diesel_power'],
name='Diesel',
marker=go.Marker(color='rgba(100, 45, 30, 1)')),
go.Bar(x=index, y=df_monthly_sums['renewable_generation_curtailed'],
name='Curtailed',
marker=go.Marker(color='rgb(220, 110, 0)'))]
captions = ""
gap = 0.12
if interval_to_average_over == 'D':
gap = 0
captions = "Daily"
if interval_to_average_over == 'W':
gap = 0.12
captions = "Weekly"
if interval_to_average_over == 'M':
gap = 0.12
captions = 'Monthly'
if interval_to_average_over == 'H':
gap = 0
captions = "Hourly"
layout = go.Layout(
barmode='stack',
bargap=gap,
plot_bgcolor='rgb(245,242,240)', # plot area color
paper_bgcolor='rgb(235,232,230)', # frame background color
height=550,
width=800,
font=go.Font(
family="Arial",
color='Black',
size=14))
layout['xaxis'].update(update_axis('<b>Annual '+ captions + ' Timeseries </b>'))
layout['yaxis'].update(update_axis('<b>Energy Mix Contribution (kWh)</b>'))
annotations = [make_anno1('<b>Energy Mix by Respective Energy Sources Summed: </b>'+ captions, 16, 0.5, 1.2)]
layout['annotations'] = annotations
fig = go.Figure(data=data, layout=layout)
if display_interactive:
po.iplot(fig)
else:
py.image.save_as(fig, format = 'jpeg',
filename = 'C:\\Users\Gregory\Desktop\Thesis Masters\RESULTS_FOLDER\PLOTS\Energy Plots\Energy_'+ \
interval_to_average_over)
In [33]:
def plot_LCOE_components(lcoe_dict):
LCOE_components = lcoe_dict
dataLCOEbar =[
go.Bar(
x = ['Battery Bank'],
y = [LCOE_components['lc_battery_cap']],
name = 'Battery Capital',
marker = go.Marker (color = 'green')),
go.Bar(
x = ['Battery Bank'],
y = [LCOE_components['lc_battery_om']],
name = 'Battery O&M',
marker = go.Marker (color = 'limegreen')),
go.Bar(
x = ['PV Array'],
y = [LCOE_components['lc_battery_cap'] + LCOE_components['lc_battery_om']],
name = '',
marker = go.Marker (color = 'rgba(1,1,1, 0.0)')),
go.Bar(
x = ['PV Array'],
y = [LCOE_components['lc_pv_cap']],
name = 'PV Capital',
marker = go.Marker (color = 'orangered')),
go.Bar(
x = ['PV Array'],
y = [LCOE_components['lc_pv_om']],
name = 'PV O&M',
marker = go.Marker (color = 'orange')),
go.Bar(
x = ['Wind Array'],
y = [LCOE_components['lc_battery_cap'] + LCOE_components['lc_battery_om'] +
LCOE_components['lc_pv_cap'] + LCOE_components['lc_pv_om']],
name = '',
marker = go.Marker (color = 'rgba(1,1,1, 0.0)')),
go.Bar(
x = ['Wind Array'],
y = [LCOE_components['lc_wind_cap']],
name = 'Wind Array Capital',
marker = go.Marker (color = 'dodgerblue')),
go.Bar(
x = ['Wind Array'],
y = [LCOE_components['lc_wind_om']],
name = 'Wind Array O&M',
marker = go.Marker (color = 'darkturquoise')),
go.Bar(
x = ['Hydrogen'],
y = [LCOE_components['lc_battery_cap'] + LCOE_components['lc_battery_om'] +
LCOE_components['lc_pv_cap'] + LCOE_components['lc_pv_om'] +
LCOE_components['lc_wind_cap'] + LCOE_components['lc_wind_om']],
name = '',
marker = go.Marker (color = 'rgba(1,1,1, 0.0)')),
go.Bar(
x = ['Hydrogen'],
y = [LCOE_components['lc_fuel_cell_cap']],
name = 'Fuel Cell Capital',
marker = go.Marker (color='rgba(255, 100, 100, 1)')),
go.Bar(
x = ['Hydrogen'],
y = [LCOE_components['lc_fuel_cell_om']],
name = 'Fuel Cell O&M',
marker = go.Marker (color='rgba(255, 120, 120, 1)')),
go.Bar(
x = ['Hydrogen'],
y = [LCOE_components['lc_electrolyzer_cap']],
name = 'Electrolyser Capital',
marker = go.Marker (color = 'rgba(0, 185, 110, 1')),
go.Bar(
x = ['Hydrogen'],
y = [LCOE_components['lc_electrolyzer_om']],
name = 'Electrolyser O&M',
marker = go.Marker (color = 'rgba(15, 200, 125, 1')),
go.Bar(
x = ['Hydrogen'],
y = [LCOE_components['lc_hydrogen_storage_cap']],
name = 'Hydrogen Storage',
marker = go.Marker (color = 'rgba(150, 200, 0, 1)')),
go.Bar(
x = ['Inverter'],
y = [LCOE_components['lc_battery_cap'] + LCOE_components['lc_battery_om'] +
LCOE_components['lc_pv_cap'] + LCOE_components['lc_pv_om'] +
LCOE_components['lc_wind_cap'] + LCOE_components['lc_wind_om'] +
LCOE_components['lc_fuel_cell_cap'] + LCOE_components['lc_fuel_cell_om'] +
LCOE_components['lc_electrolyzer_cap'] + LCOE_components['lc_electrolyzer_om'] +
LCOE_components['lc_hydrogen_storage_cap']],
name = '',
marker = go.Marker (color = 'rgba(1,1,1, 0.0)')),
go.Bar(
x = ['Inverter'],
y = [LCOE_components['lc_inverter_cap']],
name = 'Inverter Capital',
marker = go.Marker (color = 'grey')),
go.Bar(
x = ['Inverter'],
y = [LCOE_components['lc_inverter_fom']],
name = 'Inverter O&M',
marker = go.Marker (color = 'lightgrey')),
#The Below is needed until Plotly releases a native Waterfall chart
#Hack to implement it using a normal bar chart.
go.Bar(
x = ['Diesel Cost'],
y = [LCOE_components['lc_battery_cap']+LCOE_components['lc_battery_om'] +
LCOE_components['lc_pv_cap'] + LCOE_components['lc_pv_om'] +
LCOE_components['lc_wind_cap'] + LCOE_components['lc_wind_om'] +
LCOE_components['lc_fuel_cell_cap'] + LCOE_components['lc_fuel_cell_om'] +
LCOE_components['lc_electrolyzer_cap'] + LCOE_components['lc_electrolyzer_om'] +
LCOE_components['lc_hydrogen_storage_cap'] +
LCOE_components['lc_inverter_cap'] + LCOE_components['lc_inverter_fom']
],
name = '',
marker = go.Marker (color = 'rgba(1,1,1, 0.0)')),
go.Bar(
x = ['Diesel Cost'],
y = [LCOE_components['lc_diesel_cap']],
name = 'Diesel Capital',
marker = go.Marker (color = 'maroon')),
go.Bar(
x = ['Diesel Cost'],
y = [LCOE_components['lc_diesel_om']],
name = 'Diesel O&M',
marker = go.Marker (color = 'crimson')),
go.Bar(
x = ['Diesel Cost'],
y = [LCOE_components['lc_diesel_fuel']],
name = 'Diesel Fuel',
marker = go.Marker (color = 'darksalmon')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_battery_cap']],
name = 'Battery Capital',
marker = go.Marker (color = 'green')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_battery_om']],
name = 'Battery O&M',
marker = go.Marker (color = 'limegreen')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_pv_cap']],
name = 'PV Capital',
marker = go.Marker (color = 'orangered')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_pv_om']],
name = 'PV O&M',
marker = go.Marker (color = 'orange')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_wind_cap']],
name = 'Wind Capital',
marker = go.Marker (color = 'dodgerblue')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_wind_om']],
name = 'Wind O&M',
marker = go.Marker (color = 'darkturquoise')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_fuel_cell_cap']],
name = 'Fuel Cell Capital',
marker = go.Marker (color='rgba(255, 100, 100, 1)')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_fuel_cell_om']],
name = 'Fuel Cell O&M',
marker = go.Marker (color = 'rgba(255, 120, 120, 1)')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_electrolyzer_cap']],
name = 'Electrolyser Capital',
marker = go.Marker (color = 'rgba(0, 185, 110, 1')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_electrolyzer_om']],
name = 'Electrolyser O&M',
marker = go.Marker (color = 'rgba(15, 200, 125, 1')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_hydrogen_storage_cap']],
name = 'Hydrogen Storage',
marker = go.Marker (color = 'rgba(150, 200, 0, 1)')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_inverter_cap']],
name = 'Inverter Capital',
marker = go.Marker (color = 'grey')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_inverter_fom']],
name = 'Inverter O&M',
marker = go.Marker (color = 'lightgrey')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_diesel_cap']],
name = 'Diesel Capital',
marker = go.Marker (color = 'maroon')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_diesel_om']],
name = 'Diesel O&M',
marker = go.Marker (color = 'crimson')),
go.Bar(
x = ['Total Levelised Cost'],
y = [LCOE_components['lc_diesel_fuel']],
name = 'Diesel Fuel',
marker = go.Marker (color = 'darksalmon'))]
layoutLCOE = go.Layout(
barmode='stack',
bargap = 0.1,
plot_bgcolor='#EFECEA', # set plot and
paper_bgcolor='#EFECEA', # frame background color
height=550,
width=700,
font=go.Font(
family="Arial",
color='black',
size =15))
layoutLCOE['xaxis'].update(update_axis('<b>Stacked Levelised Cost Components for mini-grid</b>'))
layoutLCOE['yaxis'].update(update_axis('<b>Levelized Cost in $c/kWh</b>'))
annotationsLCOE = [make_anno1(
'<b>Levelised Costing of Mini-Grid Configuration</b><br>Component-Wise Decomposition',
18, 0.5, 1.25)]
layoutLCOE['annotations'] = annotationsLCOE
layoutLCOE['annotations'].append(
make_anno1('<b>LCOE Total:<br>'+"{:.2f}".format(run_result['LCOE'])+' $c/kWh</b>',
18, 1.45, 1.25))
figLCOE = go.Figure(data = dataLCOEbar, layout=layoutLCOE)
#po.plot(figLCOE)
po.iplot(figLCOE, show_link = False)
print("\n---Min-Grid Levelized Costing Decomposition Plotting Complete---")
In [ ]: