diff --git a/.DS_Store b/.DS_Store new file mode 100644 index 0000000..0e91fc5 Binary files /dev/null and b/.DS_Store differ diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml new file mode 100644 index 0000000..941b4f1 --- /dev/null +++ b/.github/workflows/ci.yml @@ -0,0 +1,44 @@ +# This workflow will install Python dependencies, run tests and lint with a variety of Python versions +# For more information see: https://help.github.com/actions/language-and-framework-guides/using-python-with-github-actions + +name: Python package + +on: + push: + branches: [ master ] + pull_request: + branches: [ master ] +jobs: + build: + runs-on: ubuntu-latest + strategy: + matrix: + python-version: [3.6, 3.7, 3.8, 3.9] + steps: + - uses: actions/checkout@v2 + - name: Set up Python + uses: actions/setup-python@v2 + with: + python-version: ${{ matrix.python-version }} + # You can test your matrix by printing the current Python version + - name: Display Python version + run: python -c "import sys; print(sys.version)" + - name: Install dependencies + run: | + python -m pip install --upgrade pip + pip install Cython + pip install click==8.0.4 + pip install -r requirements_dev.txt + - name: Run black + run: + black --check . + #- name: Run flake8 + # run: flake8 + #- name: Run Pylint + # run: pylint mixsimulator + - name: Run Mypy + run: mypy --allow-untyped-defs --warn-redundant-casts --strict-optional mixsimulator + #- name: tests + # run: | + # pip install .[tests] + # pytest diff --git a/.gitignore b/.gitignore old mode 100644 new mode 100755 diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml new file mode 100644 index 0000000..4bd9e29 --- /dev/null +++ b/.pre-commit-config.yaml @@ -0,0 +1,60 @@ +repos: + - repo: https://github.com/asottile/pyupgrade + rev: v2.32.1 + hooks: + - id: pyupgrade + args: [ "--py38-plus" ] + - repo: https://github.com/pre-commit/mirrors-isort + rev: v5.10.1 # Use the revision sha / tag you want to point at + hooks: + - id: isort + args: ["--profile", "black"] + - repo: https://github.com/psf/black + rev: 22.3.0 + hooks: + - id: black + #- repo: https://gitlab.com/pycqa/flake8 + # rev: 4.0.1 + # hooks: + # - id: flake8 + # language_version: python3 + - repo: https://github.com/pre-commit/pre-commit-hooks + rev: v4.2.0 + hooks: + - id: check-docstring-first + - id: check-json + - id: check-merge-conflict + - id: check-yaml + - id: debug-statements + - id: end-of-file-fixer + - id: trailing-whitespace + - id: requirements-txt-fixer + - repo: https://github.com/pre-commit/mirrors-pylint + rev: v3.0.0a4 + hooks: + - id: pylint + args: + - --max-line-length=120 + - --ignore-imports=yes + - -d duplicate-code + - repo: https://github.com/pre-commit/pygrep-hooks + rev: v1.9.0 + hooks: + - id: python-check-mock-methods + - id: python-use-type-annotations + - id: python-check-blanket-noqa + - id: python-use-type-annotations + - id: text-unicode-replacement-char + - repo: https://github.com/pre-commit/mirrors-mypy + rev: v0.950 + hooks: + - id: mypy + exclude: ^tests/ + args: + [ + #--disallow-untyped-defs, + #--check-untyped-defs, + --warn-redundant-casts, + #--no-implicit-optional, TODO! + --strict-optional + ] diff --git a/Demand.py b/Demand.py deleted file mode 100644 index b137b81..0000000 --- a/Demand.py +++ /dev/null @@ -1,99 +0,0 @@ -import pandas as pd -from prophet import Prophet -from math import pi -from math import cos, floor -import pkgutil -import csv - -class Demand: - """ - Manage the Demands data - """ - - def __init__(self, demand: float = 20, var_per_day: float = 0.1 , var_per_season: float = 0.1) -> None: - self.__var_per_day = var_per_day - self.__var_per_season = var_per_season - self.__mean_demand = demand - self.__pt = Prophet(seasonality_mode='multiplicative') - self.__periods = 12*20 - self.data_demand = None - - def set_forcast_periods(self, periods) -> None: - self.__periods = periods - - def set_data_csv(self, bind = None, raw_data = None, init_date: str = "2017-01-01", delimiter: str = ";", column: str = "Total Ventes"): - """ - The method must get a dataset with at least 3 columns - - month : int, - - year : int, - - the monthly demand in kwh (determinated by the parameter "column") - - The method also use a forcast model from prophet to predict future demand. - The periods can be set by set_forcast_periods. - """ - if raw_data is not None : - data = pd.DataFrame(raw_data) - #set columns & index - header = data.iloc[0] - data = data[1:] - data.columns = header - data = data.reset_index(drop=True) - for col in data.columns.tolist(): - try: - # convert numeric values - data[col] = pd.to_numeric(data[col]) - except: - pass - else : - data = pd.read_csv(bind, delimiter = delimiter) - - data["date"]= data ["month"].astype("str") + "-" + data ["year"].astype("str") - data["datetime"] = pd.to_datetime(data["date"]) - data_to_use = pd.DataFrame() - data_to_use[["datetime","demands"]] = data.loc[data["datetime"] >= init_date][["datetime",column]] - - #let's forcast - train = pd.DataFrame() - train[["ds","y"]] = data[["datetime",column]] - self.__pt.fit(train) - future = self.__pt.make_future_dataframe(periods = self.__periods, freq='MS') - fcst = self.__pt.predict(future) - fcst = fcst.loc[fcst["ds"] > data["datetime"].tail(1).item()] - fcst[["datetime","demands"]] = fcst[["ds","yhat"]] - self.data_demand = data_to_use.append(fcst,ignore_index=True) - - return self.data_demand - - def set_data_to(self, dataset, delimiter: str=";"): - if dataset == "Toamasina": - #by defaut we keep it "Toamasina" - data = pkgutil.get_data('mixsimulator', '/data/RIToamasina/DIR-TOAMASINA_concat.csv') - data = csv.reader(data.decode('utf-8').splitlines(), delimiter = delimiter) - self.set_data_csv(raw_data=data) - else : - #by defaut we keep it "Toamasina" - data = pkgutil.get_data('mixsimulator', '/data/RIToamasina/DIR-TOAMASINA_concat.csv') - data = csv.reader(data.decode('utf-8').splitlines(), delimiter = delimiter) - self.set_data_csv(raw_data=data) - - def get_demand(self, t): - self.data_demand.reset_index() - try : - return self.data_demand.iloc[t]["demands"]/1000 #cause data units in kwh - except IndexError : - print("Please check init_date or change the forcasting params : increase periods with Demand.set_forcast_periods") - raise - - def set_mean_demand(self, demand: float): - self.__mean_demand = demand - - def get_demand_approxima(self, t, interval): - demande = self.__mean_demand * (1 + cos(4 * pi * ( t * interval )/ 24)*self.__var_per_day + cos(2 * pi * ( t * interval ) / (24*365))* self.__var_per_season) - return demande - - def get_demand_monthly(self, t, interval): - m = t/(24*30) - m = floor(m) - # for now we divide it by 30*24 (better approximation TO DO) - demande = (self.get_demand(m)/(30*24)) * (1 + cos(4 * pi * ( t * interval )/ 24)*self.__var_per_day) - return demande diff --git a/Evaluation.py b/Evaluation.py deleted file mode 100644 index 49a6889..0000000 --- a/Evaluation.py +++ /dev/null @@ -1,340 +0,0 @@ -#import sys -#from MixSimulator import MixSimulator -import matplotlib.pyplot as plt # type: ignore -#from matplotlib.ticker import StrMethodFormatter -import numpy as np # type: ignore -from math import ceil -from math import floor -#import itertools -from typing import List -import warnings -from .nevergradBased import Optimizer as opt -from . import Demand as de -from datetime import datetime -#from matplotlib import ticker - -class EvaluationBudget: - - def __init__(self) -> None : - tmp = opt.Optimizer() - self.__available_optimizers = tmp.get_available_optimizers() - self.__marker = [',', '+', '.','<','>','p','h','H','*','x','v','^','s','1','2','3','4','8'] - - def moving_average(self, x, w): - return np.convolve(x, np.ones(w), 'valid') / w - - def set_units(self,value): - # Y units - N=value - if N>=1000000: - units=1000000 - labelY=' Millions' - else : - if N>=1000 and N<1000000: - units=1000 - labelY='k' - else : - units=1 - labelY='' - return [labelY,units] - - def plot_evaluation(self, X, Y, label_y : List['str'], label : List = ["Optimizer"], max_budgets = 0, average_wide : int = 0, plot = "save"): - #reset all possible previous plots - try : - plt.close("all") - except : - pass - - #set the moving average wide - if average_wide == 0 : - for opt_name, value in Y[label_y[0]].items(): - average_wide = ceil(len(value)/4) - #units=self.set_units(value[0]) - break - - #Label y = 1 - if len(label_y) == 1 or len(label_y) == 2 : - fig, axs = plt.subplots(2, figsize=(6, 6)) - - # data integration - it = 3 #index debut cycle - for n_axs in range(0,len(label_y)): - dict_ = Y[label_y[n_axs]] - for opt_name, value in dict_.items(): - smooth_value = self.moving_average(value,average_wide) - axs[n_axs].plot(X[(average_wide - 1):], smooth_value, marker = self.__marker[it % len(dict_)],markevery = 0.1, alpha=0.5, lw=2, label=opt_name) - #axs[0][n_axs].xaxis.set_minor_locator(ticker.MultipleLocator(len(smooth_value))) - #axs[0][n_axs].yaxis.set_major_locator(ticker.MultipleLocator(1)) - it = it + 1 - #axs[0][n_axs].xaxis.set_ticks(np.arange(min(X),max(X)+1,1.0)) - #axs[0][n_axs].yticks(np.arrange(min(smooth_value),max(smooth_value)+1,1.0)) - #label_per_opt = axs[0][n_axs].text(X[-1],smooth_value[-1],opt_name+"("+str(float("{:.2f}".format(smooth_value[-1])))+")") - #texts.append(label_per_opt) - #axs[0][n_axs].annotate( label_per_opt, - # xy = ( X[-1], smooth_value[-1]), - # xytext = (1.02*X[-1], smooth_value[-1]), - # ) - #adjust_text(texts) - - # plots parametrizations - axs[n_axs].grid() - axs[n_axs].yaxis.set_tick_params(which='major', width=1.00, length=5) - axs[n_axs].xaxis.set_tick_params(which='major', width=1.00, length=5) - axs[n_axs].xaxis.set_tick_params(which='minor', width=0.75, length=2.5, labelsize=10) - axs[n_axs].set_xlabel('Budgets') - #axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) - try : - axs[n_axs].set_ylabel(label_y[n_axs]) - except : - pass - axs[n_axs].legend() - - - for n_axs in range(0,2) : - if not axs[n_axs].has_data(): - fig.delaxes(axs[n_axs]) - - fig.tight_layout() - - if plot == "save": - try: - path = "results_"+datetime.now().strftime("%d_%m_%Y") - name = path+"/Evaluation_"+datetime.now().strftime("%d%m%Y_%H%M%S")+".png" - fig.savefig(name) - plt.show() - except FileNotFoundError: - warnings.warn("Can't find the directory "+path) - name = "Evaluation_"+datetime.now().strftime("%d%m%Y_%H%M%S")+".png" - fig.savefig(name) - plt.show() - - else : - plt.show() - - else : - #For label y more than 2 - max_col = ceil(len(label_y)/2) - min_col = floor(len(label_y)/2) - fig, axs = plt.subplots(2, max_col, figsize=(10, 8)) - - # data integration - #texts=[] - for n_axs in range(0,max_col) : - dict_ = Y[label_y[n_axs]] - it = 3 #index debut cycle - for opt_name, value in dict_.items(): - smooth_value = self.moving_average(value,average_wide) - axs[0][n_axs].plot(X[(average_wide - 1):], smooth_value, marker = self.__marker[it % len(dict_)],markevery = 0.1, alpha=0.5, lw=2, label=opt_name) - #axs[0][n_axs].xaxis.set_minor_locator(ticker.MultipleLocator(len(smooth_value))) - #axs[0][n_axs].yaxis.set_major_locator(ticker.MultipleLocator(1)) - it = it + 1 - #axs[0][n_axs].xaxis.set_ticks(np.arange(min(X),max(X)+1,1.0)) - #axs[0][n_axs].yticks(np.arrange(min(smooth_value),max(smooth_value)+1,1.0)) - #label_per_opt = axs[0][n_axs].text(X[-1],smooth_value[-1],opt_name+"("+str(float("{:.2f}".format(smooth_value[-1])))+")") - #texts.append(label_per_opt) - #axs[0][n_axs].annotate( label_per_opt, - # xy = ( X[-1], smooth_value[-1]), - # xytext = (1.02*X[-1], smooth_value[-1]), - # ) - #adjust_text(texts) - - for n_axs in range(0,min_col) : - dict_ = Y[label_y[max_col+n_axs]] - it = 3 #index debut cycle - for opt_name, value in dict_.items(): - smooth_value = self.moving_average(value,average_wide) - axs[1][n_axs].plot(X[(average_wide - 1):], smooth_value, marker = self.__marker[it % len(dict_)],markevery = 0.1, alpha=0.5, lw=2, label=opt_name) - #axs[1][n_axs].xaxis.set_minor_locator(ticker.MultipleLocator(len(smooth_value))) - #axs[1][n_axs].yaxis.set_major_locator(ticker.MultipleLocator(1)) - it = it + 1 - #axs[1][n_axs].xaxis.set_ticks(np.arange(min(X),max(X)+1,1.0)) - #axs[1][n_axs].yticks(np.arrange(min(smooth_value),max(smooth_value)+1,1.0)) - #axs[1][n_axs].annotate(opt_name+"("+str(float("{:.2f}".format(smooth_value[-1])))+")", - # xy = ( X[-1], smooth_value[-1]), - # xytext = (1.02*X[-1], smooth_value[-1]), - # ) - - - # plots parametrizations - for row in range (0,2): - if row == 0 : - for n_axs in range(0,max_col) : - axs[row][n_axs].grid() - axs[row][n_axs].yaxis.set_tick_params(which='major', width=1.00, length=5) - axs[row][n_axs].xaxis.set_tick_params(which='major', width=1.00, length=5) - axs[row][n_axs].xaxis.set_tick_params(which='minor', width=0.75, length=2.5, labelsize=10) - axs[row][n_axs].set_xlabel('Budgets') - #axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) - axs[row][n_axs].set_ylabel(label_y[n_axs]) - axs[row][n_axs].legend() - else : - for n_axs in range(0,min_col) : - axs[row][n_axs].grid() - axs[row][n_axs].yaxis.set_tick_params(which='major', width=1.00, length=5) - axs[row][n_axs].xaxis.set_tick_params(which='major', width=1.00, length=5) - axs[row][n_axs].xaxis.set_tick_params(which='minor', width=0.75, length=2.5, labelsize=10) - axs[row][n_axs].set_xlabel('Budgets') - #axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) - axs[row][n_axs].set_ylabel(label_y[max_col+n_axs]) - axs[row][n_axs].legend() - - for row in range (0,2): - for n_axs in range(0,max_col) : - if not axs[row][n_axs].has_data(): - fig.delaxes(axs[row][n_axs]) - - fig.tight_layout() - - if plot == "save": - try: - path = "results_"+datetime.now().strftime("%d_%m_%Y") - name = path+"/Evaluation_"+datetime.now().strftime("%d%m%Y_%H%M%S")+".png" - fig.savefig(name) - plt.show() - except FileNotFoundError: - warnings.warn("Can't find the directory "+path) - name = "Evaluation_"+datetime.now().strftime("%d%m%Y_%H%M%S")+".png" - fig.savefig(name) - plt.show() - - else : - plt.show() - - def plot_time_evolution(self, data, label_y : List['str'], label : List = ["Optimizer"], max_budgets = 0): - #init subplot - print(data) - - fig, axs = plt.subplots(1, len(label_y)-1, figsize=(12, 4)) - - # data integration - for n_axs in range(1,len(axs)) : - dict_ = data[label_y[n_axs]] - for opt_name, value in dict_.items(): - axs[n_axs].scatter(data["execution_time (s)"][opt_name], value, alpha=0.5, lw=2, label=str(opt_name)) - - # plots parametrizations - for n_axs in range(0,len(axs)) : - axs[n_axs].grid() - axs[n_axs].yaxis.set_tick_params(length=0) - axs[n_axs].xaxis.set_tick_params(length=0) - axs[n_axs].set_xlabel('Budgets') - #axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) - axs[n_axs].set_ylabel(label_y[n_axs]) - axs[n_axs].legend() - - fig.tight_layout() - plt.show() - - def check_opt_list(self,optimizer_list): - for index in range(0, len(optimizer_list)): - try : - if optimizer_list[index] not in self.__available_optimizers: - optimizer_list.pop(index) - except IndexError : - self.check_opt_list(optimizer_list) - - - def evaluate(self, mix, sequence, max_budgets, optimizer_list: List['str'], indicator_list: List['str'], num_worker: int = 1, bind: str = None, time_index: int = 24, carbonProdLimit: float = 39500000000, time_interval : float = 1, average_wide : int = 0, penalisation : float = 1000000000000, carbon_cost: float = 0) : - #setting dataset - if bind != None: - mix.set_data_csv(str(bind)) - - self.check_opt_list(optimizer_list) - if optimizer_list == [] : - raise IndexError("Selected optimizers are not available.") - - #process - y_tmp = {} - budget = np.arange(0, max_budgets, sequence) - - ind_per_opt = {} - for opt_name in optimizer_list: - opt_index = opt.Optimizer(opt=[opt_name], budget = [max_budgets], num_worker = num_worker) - data = mix.optimizeMix(carbon_quota = carbonProdLimit, - time_interval = time_interval, optimizer = opt_index, step = sequence, penalisation = penalisation, time_index = time_index, carbon_cost = carbon_cost, plot = "save") - ind_per_opt.update({opt_name:data}) - - for indicator in indicator_list: - new_ind_per_opt = {} - for opt_name, values in ind_per_opt.items(): - ind_per_budget = [] - for budget_value in values: - ind_per_budget.append(float(budget_value[indicator])) - new_ind_per_opt.update({opt_name:ind_per_budget}) - y_tmp.update({indicator: new_ind_per_opt}) - - #plotting - self.plot_evaluation(X=np.array(budget),Y=y_tmp,label_y = indicator_list, label=optimizer_list, max_budgets = max_budgets, average_wide = average_wide) - - #return X, Y, opt_list, max_budgets - return [np.array(budget),y_tmp,optimizer_list,max_budgets] - - # NOT USE : NEED VERIFICATION - def evaluate_total_time(self, mix, sequence, max_budgets, optimizer_list: List['str'], - indicator_list: List['str'], bind = None, carbonProdLimit: float = 500000, - time_index: int = 24*265, time_interval : float = 1, average_wide : int = 0, penalisation : float = 1000000000000, plot : str = "default"): - #setting dataset - - budget = np.arange(0, max_budgets, sequence) - - if bind != None: - mix.set_data_csv(str(bind)) - - self.check_opt_list(optimizer_list) - if optimizer_list == [] : - raise IndexError("Selected optimizers are not available.") - - #process - data_interval = [] - current_demand=de.Demand(mix.get_demand(),0.2,0.2) - for time in range(0,time_index): - mix.set_demand(current_demand.get_demand_approxima(time,time_interval)) - ind_per_opt = {} - for opt_name in optimizer_list: - data = mix.optimizeMix(carbonProdLimit= carbonProdLimit, - time_interval = time_interval, optimize_with = [opt_name], budgets = [max_budgets], step = sequence, penalisation = penalisation) - ind_per_opt.update({opt_name:data}) - data_interval.append(ind_per_opt) - - Y = [] - indicator_list_WO_penalisation = indicator_list.copy() - try: - indicator_list_WO_penalisation.remove("penalized production cost (loss)") - except: - pass - for time in range(0,time_index): - y_tmp = {} - for indicator in indicator_list_WO_penalisation: - new_ind_per_opt = {} - for opt_name, values in ind_per_opt.items(): - ind_per_budget = [] - for budget_value in values: - ind_per_budget.append(float(budget_value[indicator])) - new_ind_per_opt.update({opt_name:ind_per_budget}) - y_tmp.update({indicator: new_ind_per_opt}) - Y.append(y_tmp) - - result = {} - for indicator in indicator_list: - optimizers_dict = {} - for opt_name in optimizer_list: - per_budget = [] - for budget_step in range(0,len(budget)): - value = 0. - for time in range (0, time_index): - if indicator == "penalized production cost (loss)": - value += Y[time]["production_cost ($)"][opt_name][budget_step] + np.abs((mix.get_penalisation_cost() * Y[time]["unsatisfied_demand (MWh)"][opt_name][budget_step])) - else : - value += Y[time][indicator][opt_name][budget_step] - if indicator == "penalized production cost (loss)": - value = np.log10(value) - per_budget.append(value) - optimizers_dict.update({opt_name:per_budget}) - result.update({indicator: optimizers_dict}) - - #plotting - self.plot_evaluation(X=np.array(budget),Y=result,label_y = indicator_list, label=optimizer_list, max_budgets = max_budgets,average_wide = average_wide, plot = plot) - - #return X, Y, opt_list, max_budgets - return [np.array(budget),result,optimizer_list,max_budgets] - \ No newline at end of file diff --git a/LICENSE b/LICENSE old mode 100644 new mode 100755 diff --git a/MANIFEST.in b/MANIFEST.in old mode 100644 new mode 100755 diff --git a/MixSimulator.py b/MixSimulator.py deleted file mode 100644 index f0580aa..0000000 --- a/MixSimulator.py +++ /dev/null @@ -1,382 +0,0 @@ -from .centrals import PowerCentral as pc -from .centrals import HydroCentral as hc -from .Demand import Demand -import nevergrad as ng -from .nevergradBased.Optimizer import Optimizer -import numpy as np # type: ignore -import pandas as pd # type: ignore -import pkgutil -import csv -import os -import warnings -from math import ceil -#import time -from typing import List, Any, Type, Dict -from datetime import datetime -import matplotlib.pyplot as plt # type: ignore - -class MixSimulator: - """ - The simulator Base - """ - def __init__(self, carbon_cost: float = 0, penalisation_cost: float = 1000000000000) -> None: - self.__centrals : List[Any] = [] - self.__reset_centrals() - self.__demand = Demand(20, 0.2, 0.2) - self.__lost = 0. - self.__penalisation_cost = penalisation_cost - self.__optimizer = Optimizer() - self.__carbon_cost = carbon_cost - self.__carbon_quota = 800139. # (g/MWh) - - # for reuse and get_params() - self.time_index = 24*7 - self.step = 1 - self.time_interval = 1 - self.plot = "default" - self.average_wide = 0 - - def __reset_centrals(self): - self.__centrals = [] - - ## DATA ## - def set_variation_csv(self, bind = None, delimiter: str=";") -> None: - if self.__centrals == []: - warnings.warn("Please load a original dataset by using MixSimulator.set_data_csv(...)") - raise - else : - try : - data = pd.DataFrame(pd.read_csv(bind,delimiter=delimiter)) - except FileNotFoundError as e : - print("Error occured on pandas.read_csv : ",e) - print("Please check your file") - raise - except Exception as e: - print("Error occured on pandas.read_csv : ",e) - raise - - for central_index in range(len(self.__centrals)): - if self.__centrals[central_index].is_tuneable(): - for i in range (0,data.shape[0]): - if self.__centrals[central_index].get_id() == data["centrals"][i]: - tmp_list=[] - upper = str(data["upper"][i]).split(":") - upper = [float(numeric_string) for numeric_string in upper] - lower = str(data["lower"][i]).split(":") - lower = [float(numeric_string) for numeric_string in lower] - discrete = str(data["discrete"][i]).split(":") - discrete = [float(numeric_string) for numeric_string in discrete] - self.__centrals[central_index].set_variation_params(lower = lower, upper = upper, choices = discrete) - - def set_data_csv(self, bind = None, raw_data = None, delimiter: str=";"): - if raw_data is not None : - data = pd.DataFrame(raw_data) - #set columns & index - header = data.iloc[0] - data = data[1:] - data.columns = header - data = data.reset_index(drop=True) - for column in data.columns.tolist(): - try: - # convert numeric values - data[column] = pd.to_numeric(data[column]) - except: - pass - else : - try : - data = pd.DataFrame(pd.read_csv(bind,delimiter=delimiter)) - except FileNotFoundError as e : - print("Error occured on pandas.read_csv : ",e) - print("Please check your file") - raise - except Exception as e: - print("Error occured on pandas.read_csv : ",e) - raise - - self.__reset_centrals() - centrale = pc.PowerCentral() - try : - for i in range (0,data.shape[0]): - isHydro = data["hydro"][i] - if isHydro == True : - centrale = hc.HydroCentral(data["height"][i],data["flow"][i],data["capacity"][i],data["stock_available"][i],0.1,0.8) - else : - centrale = pc.PowerCentral() - centrale.set_tuneable(data["tuneable"][i]) - centrale.set_id(str(data["centrals"][i])) - centrale.set_fuel_consumption(data["fuel_consumption"][i]) - centrale.set_availability(data["availability"][i]) - centrale.set_fuel_cost(data["fuel_cost"][i]) - centrale.set_initial_value(data["init_value"][i]) - centrale.set_lifetime(data["lifetime"][i]) - centrale.set_carbon_prod(data["carbon_production"][i]) - centrale.set_raw_power(data["raw_power"][i]) - centrale.set_nb_employees(data["nb_employees"][i]) - centrale.set_mean_employees_salary(data["mean_salary"][i]) - centrale.set_max_var(data["max_var"][i]) - self.__centrals.append(centrale) - self.__demand.set_mean_demand(data["Demand"][0]) - self.__lost=data["lost"][0] - except KeyError: - print("Columns must be in: tuneable, centrals, fuel_consumption, availability, fuel_cost, init_value, lifetime, carbon_cost, raw_power, nb_employees, mean_salary, demand, lost, height, flow, capacity, stock_available") - raise - - def set_data_to(self, dataset, delimiter: str=";"): - if dataset == "Toamasina": - #by defaut we keep it "Toamasina" - data = pkgutil.get_data('mixsimulator', '/data/RIToamasina/dataset_RI_Toamasina_v2.csv') - data = csv.reader(data.decode('utf-8').splitlines(), delimiter = delimiter) - self.set_data_csv(raw_data=data) - else : - #by defaut we keep it "Toamasina" - data = pkgutil.get_data('mixsimulator', '/data/RIToamasina/dataset_RI_Toamasina_v2.csv') - data = csv.reader(data.decode('utf-8').splitlines(), delimiter = delimiter) - self.set_data_csv(raw_data=data) - - def set_demand(self, demand: Demand): - self.__demand = demand - - def get_demand(self): - return self.__demand - - def set_lost(self, lost: float): - self.__lost = lost - - def set_optimizer(self, optimizer: Optimizer): - self.__optimizer = optimizer - - def get_optimizer(self) -> Optimizer: - return self.__optimizer - - def set_carbon_quota(self, cb_quota: float ): - self.__carbon_quota = cb_quota - - def set_penalisation_cost(self, k): - self.__penalisation_cost = k - - def get_penalisation_cost(self): - return self.__penalisation_cost - - def set_carbon_cost(self, carbon_cost): - self.__carbon_cost = carbon_cost - - def get_carbon_cost(self): - return self.__carbon_cost - - ## EVALUATION FONCTION ## - def get_production_cost_at_t(self, usage_coef, time_index, time_interval): - production_cost = 0 - for centrale_index in range (0, len(self.__centrals)): - central = self.__centrals[centrale_index] - production_cost += ( (central.get_fuel_cost() * central.get_fuel_consumption() - * central.get_raw_power() * usage_coef[centrale_index]) + ( central.get_employees_salary() - + central.get_amortized_cost(time_index) ) ) * time_interval - return production_cost - - def get_production_at_t(self, usage_coef, time_interval): - production = 0 - for centrale_index in range (0, len(self.__centrals)): - central = self.__centrals[centrale_index] - production += (central.get_raw_power() * usage_coef[centrale_index]) * time_interval - return production - - def get_unsatisfied_demand_at_t(self, usage_coef, time_index, time_interval): - #return ( self.__demand.get_demand_approxima(time_index, time_interval) - self.get_production_at_t(usage_coef, time_interval)) - return ( self.__demand.get_demand_monthly(time_index, time_interval) - self.get_production_at_t(usage_coef, time_interval)) - - def get_carbon_production_at_t(self, usage_coef, time_interval): - carbon_prod = 0 - for centrale_index in range (0, len(self.__centrals)): - central = self.__centrals[centrale_index] - carbon_prod += (central.get_carbon_production() * usage_coef[centrale_index] * central.get_raw_power()) * time_interval - return carbon_prod - - def get_carbon_over_production(self, usage_coef, time_interval): - emited_carbon = 0 # (g) - total_production = 1 - for t in range(0, len(usage_coef)): - emited_carbon += self.get_carbon_production_at_t(usage_coef[t], time_interval) - total_production += self.get_production_at_t(usage_coef[t], time_interval) - carbon_production = emited_carbon/total_production - return max(0, carbon_production - self.__carbon_quota) # (g/MWh) - - def get_weighted_coef(self, usage_coef, time_interval): - for central_index in range(0, len(usage_coef[0])): - self.__centrals[central_index].reset_central() - weighted_coef = usage_coef.copy() - for t in range(0, len(weighted_coef)): - for central_index in range(0, len(weighted_coef[t])): - min_av = self.__centrals[central_index].get_min_availability(t) - max_av = self.__centrals[central_index].get_max_availability(t) - if max_av < min_av: - min_av = 0 # a verifier - weighted_coef[t][central_index] = min_av + weighted_coef[t][central_index]*(max_av-min_av) - self.__centrals[central_index].back_propagate(weighted_coef[t][central_index], t, time_interval) - return weighted_coef - - def loss_function(self, usage_coef, time_interval : int = 1) -> float : - usage_coef = self.arrange_coef_as_array_of_array(usage_coef) - weighted_coef = self.get_weighted_coef(usage_coef, time_interval=time_interval) - loss = 0 - for t in range(0, len(weighted_coef)): - loss += self.get_production_cost_at_t(weighted_coef[t], t, time_interval) + ( self.get_penalisation_cost() * np.abs( self.get_unsatisfied_demand_at_t(weighted_coef[t], t, time_interval)) ) - loss += self.get_carbon_cost() * (self.get_carbon_over_production(weighted_coef, time_interval) ) - return loss - - def arrange_coef_as_array_of_array(self, raw_usage_coef): - ordered_coef = [] - cur_time_coef = [] - for coef_index in range(len(raw_usage_coef)): - cur_time_coef.append(raw_usage_coef[coef_index]) - ## indice de la premiere centrale a t+1 - if (coef_index+1)%len(self.__centrals) == 0: - ordered_coef.append(cur_time_coef) - cur_time_coef = [] - return ordered_coef - - - ## OPTiMiZATION ## - - def __opt_params(self, time_index): - self.__optimizer.set_parametrization(self.get_opt_params(time_index)) - - def get_opt_params(self, time_index): - variable_parametrization = [] - for _ in range(time_index): - for central_index in range(len(self.__centrals)): - if not self.__centrals[central_index].is_tuneable(): - variable_parametrization += [ng.p.Choice([0.,1.])] - else: - #check the params by - #print(self.__centrals[central_index].get_variation_params()) - variable_parametrization += [self.__centrals[central_index].get_variation_params()] - return ng.p.Tuple(*variable_parametrization) - - def optimizeMix(self, carbon_quota: float = None, demand: Demand = None, lost: float = None, - optimizer: Optimizer = None, step : int = 1, - time_index: int = 24*7, time_interval: float = 1, - penalisation : float = None, carbon_cost : float = None, plot : str = "default", average_wide : int = 0): - - # init params - self.time_index = time_index - self.step = step - self.time_interval = time_interval - self.plot = plot - self.average_wide = average_wide - if demand is not None : self.set_demand(demand) - if lost is not None : self.set_lost(lost) - if penalisation is not None : self.set_penalisation_cost(penalisation) - if carbon_cost is not None : self.set_carbon_cost(carbon_cost) - if carbon_quota is not None : self.set_carbon_quota(carbon_quota) - if optimizer is not None : self.set_optimizer(optimizer) - - # tune optimizer parametrization - self.__opt_params(self.time_index) - - #init constraints - constraints = {} - constraints.update({"time_interval":self.time_interval}) - - #let's optimize - results = self.__optimizer.optimize(mix = self , func_to_optimize = self.loss_function, constraints=constraints, step = self.step) - - results = self.__reshape_results(results, self.time_interval) - - self.plotResults(results, mode = self.plot , time_interval = self.time_interval, average_wide = self.average_wide) - - return results - - def __reshape_results(self, results, time_interval): - for tmp in results: - usage_coef = self.arrange_coef_as_array_of_array(tmp['coef']) - tmp.update({"coef":self.get_weighted_coef(usage_coef, time_interval)}) - - # for central_index in range(0, len(self.__centrals)): - # try: - # print(self.__centrals[central_index].get_stock_evolution()) - # # Not a hydro power plant, so the methode does not exist - # except: - # pass - return results - - def get_params(self) -> dict: - return {"centrals" : self.__centrals, "optimizer" : self.get_optimizer(), - "penalisation_cost" : self.get_penalisation_cost(), "carbon_cost" : self.get_carbon_cost(), - "demand" : self.__demand, "lost" : self.__lost, "carbon_quota" : self.__carbon_quota, - "step" : self.step, "time_interval" : self.time_interval, "time_index" : self.time_index, - "plot" : self.plot, "moving average_wide" : self.average_wide} - - ## PLOT ## - def moving_average(self, x, w): - return np.convolve(x, np.ones(w), 'valid') / w - - def plotResults(self, optimum : dict = {} , mode : str = "default", time_interval : float = 1, average_wide : int = 0): - #set the moving average wide - if average_wide == 0 : - average_wide = ceil(len(optimum[-1]["coef"])/4) - - if mode == "default" or mode == "save": - #set Y - Y: Dict[str,List[float]] ={} - label_y: List[str]=[] - for c in self.__centrals : - label_y.append(c.get_id()) - Y.update({c.get_id():[]}) - for array in optimum[-1]["coef"] : - for enum, value in enumerate(array) : - Y[label_y[enum]].append(value) - - fig, axs = plt.subplots(1, 1, figsize=(6, 6)) - - # data integration - X = [i for i in range(0,len(optimum[-1]["coef"]))] - for n_axs in range(0,1) : - for central, value in Y.items(): - smooth_value = self.moving_average(value,average_wide) - axs.plot(X[(average_wide - 1):], smooth_value, '.-' ,alpha=0.5, lw=2, label=central) - - - # Add execution_time and loss information - info = "production_cost: "+ "{:.3f}".format(optimum[-1]["loss"])+" ; execution_time: "+"{:.3f}".format(optimum[-1]["elapsed_time"]) - plt.annotate(info, - xy=(0.5, 0), xytext=(0, 10), - xycoords=('axes fraction', 'figure fraction'), - textcoords='offset points', - size=10, ha='center', va='bottom') - - # plots parametrizations - axs.grid() - axs.yaxis.set_tick_params(length=0) - axs.xaxis.set_tick_params(length=0) - axs.set_xlabel('hours') - #axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) - axs.set_ylabel('coef (%)') - axs.legend() - - fig.tight_layout() - try : - path = "results_"+datetime.now().strftime("%d_%m_%Y") - os.makedirs(path) - name = path+"/"+"opt_"+str(self.get_optimizer().get_optimizers())+"_"+datetime.now().strftime("%H%M%S")+".png" - fig.savefig(name) - if mode == "default" : - plt.show() - - except OSError: - warnings.warn("Can't create the directory "+path+" or already exists") - try : - name = path+"/"+"opt_"+str(self.get_optimizer().get_optimizers())+"_"+datetime.now().strftime("%H%M%S")+".png" - fig.savefig(name) - if mode == "default" : - plt.show() - except FileNotFoundError: - name = "opt_"+str(self.get_optimizer().get_optimizers())+"_"+datetime.now().strftime("%H%M%S")+".png" - fig.savefig(name) - if mode == "default" : - plt.show() - - elif mode == "None" : - pass - else : - warnings.warn("Choose an available option : default, save or None") - #plt.show() diff --git a/README.md b/README.md old mode 100644 new mode 100755 index ea7d0dc..3733ccb --- a/README.md +++ b/README.md @@ -1,15 +1,22 @@ -# MixSimulator +# MixSimulator MixSimulator is an application with an optimization model for calculating and simulating the least cost of an energy mix under certain constraints. The optimizers used are based on the Nevergrad Python 3.6+ library. The primary objective of the simulator is to study the relevance of an energy mix connected to each Inter-connected Grid through the coefficient of usage of each unit in the production cost. +## Version 0.4 +The current version is a multi-agent system (MAS) approach but keeps the previous classic optimization approach available. Check `test_mas.py` for more details. (Available events are : powerplant shutdown and powerplant power_on). + +An experiment on evaluating both method is available in `Experiments/Scenario_type.py` or [![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/drive/1IzjJlNg3fCp14GICB2GGYYwmSIEK9yQf?usp=sharing) + +#### Note +This project is a work in progress so it can not yet used in production (Many changes are on their way). Feedbacks are welcome! + ### Specifications : - Generic simulator, compatible with data from Madagascar and those from abroad (but may require data pre-processing beforehand); - Optimization of duty cycle (or usage coefficient) under constraints ; - Get the production cost and various performance indicators (CO2 emission, unsatisfied demand); - Estimate of the costs of a mix or a power plant over the long term ; -- Comparison of the performance indicators on different optimizers. - +- Comparison of the performance indicators on different optimizers. (see `EvaluationBudget` method) ### Perspectives : - Add other constraints ; @@ -26,15 +33,7 @@ pip install mixsimulator MixSimulator is written in Python 3.6 and requires the following Python packages : nevergrad, prophet, typing, numpy, pandas and matplotlib. ## How to run -As MixSimulator is a python package, it can be called and used as we can see in `test.py`. - -List of classes and directories : -- MixSimulator : System basis (Adaptation of the Nevergrad optimizers to the project and auto-parameterization) ; -- nevergradBased/Optimizer : Instance of the optimizer (Receives the indications on the optimizer to choose and the parameters); -- centrals/* : Gathers all the common specifications of the control units (central) ; -- Evaluation : Class for evaluating mix based on performance indicators on several optimizers ; -- data/ : Groups the available datasets. -- documentation/ : documents about the project. +As MixSimulator is a python package, it can be called and used as we can see in `test_classic.py` and `test_mas.py`. Official documentation will accompany the first release version. @@ -70,13 +69,22 @@ Hydro specification : `nb_employees * mean_salary` **can be used as a variable cost of the plant if you want to directly use other informations as variable cost.** ### Demand and Variation datas -There is also "DIR-TOAMASINA_concat.csv" about Consumption data (in kwh, more details in Demand.py) and "dataset_RI_Toamasina_variation_template.csv" about limits in variation of power plants load following (WIP). +There is also "DIR-TOAMASINA_concat.csv" about Consumption data (in kwh, more details in `demand/`) and "dataset_RI_Toamasina_variation_template.csv" about limits in variation of power plants load following (WIP). **If you have datasets of any region in the world that can be used to evaluate our model, please contact us.** +## Citation +``` +@misc{MixSimulator, + author = {Solofohanitra R. Andriamalala and Toky A. Andriamizakason and Andry Rasoanaivo}, + title = {{MixSimulator - An electricity mix simulator and optimization}}, + year = {2020}, + publisher = {GitHub}, + journal = {GitHub repository}, + howpublished = {\url{https://GitHub.com/Foloso/MixSimulator}}, +} +``` +You can also read the associated paper : [Multi-agent vs Classic System of an Electricity Mix Production Optimization](https://link.springer.com/chapter/10.1007/978-3-031-30229-9_8) ## Contact For questions and feedbacks related to the project, please send an email to r.andry.rasoanaivo@gmail.com or soloforahamefy@gmail.com or tokyandriaxel@gmail.com - -## Note -This project is a work in progress so it can not yet used in production (Many changes are on their way). Feedbacks are welcome! diff --git a/centrals/__init__.pyc b/centrals/__init__.pyc deleted file mode 100644 index febd5b9..0000000 Binary files a/centrals/__init__.pyc and /dev/null differ diff --git a/data/RIToamasina/dataset_RI_Toamasina_variation_template.csv b/data/RIToamasina/dataset_RI_Toamasina_variation_template.csv deleted file mode 100644 index 6149331..0000000 --- a/data/RIToamasina/dataset_RI_Toamasina_variation_template.csv +++ /dev/null @@ -1,4 +0,0 @@ -centrals;lower;upper;discrete -572;0.5;1;0:0.25 -573;0;0.5;0.75:1 -Betainomby;0;0;0:0.25:0.5:0.75:1 \ No newline at end of file diff --git a/demand_simulator/demand_simulator.py b/demand_simulator/demand_simulator.py deleted file mode 100644 index 2d1a138..0000000 --- a/demand_simulator/demand_simulator.py +++ /dev/null @@ -1,28 +0,0 @@ -import pandas as pd -import numpy as np -from sklearn.linear_model import LinearRegression - - -def get_accuracy(original_score, pred): - return ((original_score-pred)**2).sum() - - -data = pd.read_csv("../data/RIToamasina/demand/toamasina_demand.csv") -columns = ['MT/HT Abonnes', 'BT Abonnes', 'Total Abonnes', 'year', 'month'] -train = data.loc[data["year"]!=2017] -test = data.loc[data["year"]==2017] - -X_train = train[columns] -y_1_train = train["Production Brute Totale"] -y_2_train = train["Total Ventes"] - -X_test = test[columns] -y_1_test = test["Production Brute Totale"] -y_2_test = test["Total Ventes"] - - -reg = LinearRegression().fit(X_train, y_1_train) - -pred = reg.predict(X_test) - -print(get_accuracy(pred, y_1_test)) \ No newline at end of file diff --git a/documentation/Project_Documentation.pdf b/documentation/Project_Documentation.pdf old mode 100644 new mode 100755 diff --git a/documentation/images/Evolution_coef_72h.png b/documentation/images/Evolution_coef_72h.png old mode 100644 new mode 100755 diff --git a/mixsimulator/Evaluation.py b/mixsimulator/Evaluation.py new file mode 100644 index 0000000..edaded6 --- /dev/null +++ b/mixsimulator/Evaluation.py @@ -0,0 +1,300 @@ +import warnings +from datetime import datetime +from math import ceil, floor +from typing import Any, List + +import matplotlib.pyplot as plt # type: ignore + +# from matplotlib.ticker import StrMethodFormatter +import numpy as np # type: ignore + +from .demand.classic import Demand as de +from .nevergradBased import Optimizer as opt + + +class EvaluationBudget: + def __init__(self) -> None: + tmp = opt.Optimizer() + self.__available_optimizers = tmp.get_available_optimizers() + self.__marker = [",", "+", ".", "<", ">", "p", "h", "H", "*", "x", "v", "^", "s", "1", "2", "3", "4", "8"] + + def moving_average(self, x, w): + return np.convolve(x, np.ones(w), "valid") / w + + def set_units(self, value): + # Y units + N = value + if N >= 1000000: + units = 1000000 + labelY = " Millions" + else: + if N >= 1000 and N < 1000000: + units = 1000 + labelY = "k" + else: + units = 1 + labelY = "" + return [labelY, units] + + def plot_evaluation( + self, + X, + Y, + label_y: List["str"], + label: List = ["Optimizer"], + max_budgets=0, + average_wide: int = 0, + plot="save", + ): + # reset all possible previous plots + try: + plt.close("all") + except: + pass + + # set the moving average wide + if average_wide == 0: + for opt_name, value in Y[label_y[0]].items(): + average_wide = ceil(len(value) / 4) + # units=self.set_units(value[0]) + break + + # Label y = 1 + if len(label_y) == 1 or len(label_y) == 2: + fig, axs = plt.subplots(2, figsize=(6, 6)) + + # data integration + it = 3 # index debut cycle + for n_axs in range(0, len(label_y)): + dict_ = Y[label_y[n_axs]] + for opt_name, value in dict_.items(): + smooth_value = self.moving_average(value, average_wide) + axs[n_axs].plot( + X[(average_wide - 1) :], + smooth_value, + marker=self.__marker[it % len(dict_)], + markevery=0.1, + alpha=0.5, + lw=2, + label=opt_name, + ) + it = it + 1 + + # plots parametrizations + axs[n_axs].grid() + axs[n_axs].yaxis.set_tick_params(which="major", width=1.00, length=5) + axs[n_axs].xaxis.set_tick_params(which="major", width=1.00, length=5) + axs[n_axs].xaxis.set_tick_params(which="minor", width=0.75, length=2.5, labelsize=10) + axs[n_axs].set_xlabel("Budgets") + # axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) + try: + axs[n_axs].set_ylabel(label_y[n_axs]) + except: + pass + axs[n_axs].legend() + + for n_axs in range(0, 2): + if not axs[n_axs].has_data(): + fig.delaxes(axs[n_axs]) + + fig.tight_layout() + + if plot == "save": + try: + path = "results_" + datetime.now().strftime("%d_%m_%Y") + name = path + "/Evaluation_" + datetime.now().strftime("%d%m%Y_%H%M%S") + ".png" + fig.savefig(name) + plt.show() + except FileNotFoundError: + warnings.warn("Can't find the directory " + path) + name = "Evaluation_" + datetime.now().strftime("%d%m%Y_%H%M%S") + ".png" + fig.savefig(name) + plt.show() + + else: + plt.show() + + else: + # For label y more than 2 + max_col = ceil(len(label_y) / 2) + min_col = floor(len(label_y) / 2) + fig, axs = plt.subplots(2, max_col, figsize=(10, 8)) + + # data integration + # texts=[] + for n_axs in range(0, max_col): + dict_ = Y[label_y[n_axs]] + it = 3 # index debut cycle + for opt_name, value in dict_.items(): + smooth_value = self.moving_average(value, average_wide) + axs[0][n_axs].plot( + X[(average_wide - 1) :], + smooth_value, + marker=self.__marker[it % len(dict_)], + markevery=0.1, + alpha=0.5, + lw=2, + label=opt_name, + ) + it = it + 1 + + for n_axs in range(0, min_col): + dict_ = Y[label_y[max_col + n_axs]] + it = 3 # index debut cycle + for opt_name, value in dict_.items(): + smooth_value = self.moving_average(value, average_wide) + axs[1][n_axs].plot( + X[(average_wide - 1) :], + smooth_value, + marker=self.__marker[it % len(dict_)], + markevery=0.1, + alpha=0.5, + lw=2, + label=opt_name, + ) + it = it + 1 + + # plots parametrizations + for row in range(0, 2): + if row == 0: + for n_axs in range(0, max_col): + axs[row][n_axs].grid() + axs[row][n_axs].yaxis.set_tick_params(which="major", width=1.00, length=5) + axs[row][n_axs].xaxis.set_tick_params(which="major", width=1.00, length=5) + axs[row][n_axs].xaxis.set_tick_params(which="minor", width=0.75, length=2.5, labelsize=10) + axs[row][n_axs].set_xlabel("Budgets") + # axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) + axs[row][n_axs].set_ylabel(label_y[n_axs]) + axs[row][n_axs].legend() + else: + for n_axs in range(0, min_col): + axs[row][n_axs].grid() + axs[row][n_axs].yaxis.set_tick_params(which="major", width=1.00, length=5) + axs[row][n_axs].xaxis.set_tick_params(which="major", width=1.00, length=5) + axs[row][n_axs].xaxis.set_tick_params(which="minor", width=0.75, length=2.5, labelsize=10) + axs[row][n_axs].set_xlabel("Budgets") + # axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) + axs[row][n_axs].set_ylabel(label_y[max_col + n_axs]) + axs[row][n_axs].legend() + + for row in range(0, 2): + for n_axs in range(0, max_col): + if not axs[row][n_axs].has_data(): + fig.delaxes(axs[row][n_axs]) + + fig.tight_layout() + + if plot == "save": + try: + path = "results_" + datetime.now().strftime("%d_%m_%Y") + name = path + "/Evaluation_" + datetime.now().strftime("%d%m%Y_%H%M%S") + ".png" + fig.savefig(name) + plt.show() + except FileNotFoundError: + warnings.warn("Can't find the directory " + path) + name = "Evaluation_" + datetime.now().strftime("%d%m%Y_%H%M%S") + ".png" + fig.savefig(name) + plt.show() + + else: + plt.show() + + def plot_time_evolution(self, data, label_y: List["str"], label: List = ["Optimizer"], max_budgets=0): + # init subplot + print(data) + + fig, axs = plt.subplots(1, len(label_y) - 1, figsize=(12, 4)) + + # data integration + for n_axs in range(1, len(axs)): + dict_ = data[label_y[n_axs]] + for opt_name, value in dict_.items(): + axs[n_axs].scatter(data["execution_time (s)"][opt_name], value, alpha=0.5, lw=2, label=str(opt_name)) + + # plots parametrizations + for n_axs in range(0, len(axs)): + axs[n_axs].grid() + axs[n_axs].yaxis.set_tick_params(length=0) + axs[n_axs].xaxis.set_tick_params(length=0) + axs[n_axs].set_xlabel("Budgets") + # axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) + axs[n_axs].set_ylabel(label_y[n_axs]) + axs[n_axs].legend() + + fig.tight_layout() + plt.show() + + def check_opt_list(self, optimizer_list): + for index in range(0, len(optimizer_list)): + try: + if optimizer_list[index] not in self.__available_optimizers: + optimizer_list.pop(index) + except IndexError: + self.check_opt_list(optimizer_list) + + def evaluate( + self, + mix, + sequence, + max_budgets, + optimizer_list: List["str"], + indicator_list: List["str"], + num_worker: int = 1, + bind: Any = None, + time_index: int = 24, + carbonProdLimit: float = 39500000000, + time_interval: float = 1, + average_wide: int = 0, + penalisation: float = 1000000000000, + carbon_cost: float = 0, + ): + # setting dataset + if bind: + mix.set_data_csv(str(bind)) + else: + mix.set_data_to("Toamasina") + + self.check_opt_list(optimizer_list) + if optimizer_list == []: + raise IndexError("Selected optimizers are not available.") + + # process + y_tmp = {} + budget = np.arange(0, max_budgets, sequence) + + ind_per_opt = {} + for opt_name in optimizer_list: + opt_index = opt.Optimizer(opt=[opt_name], budget=[max_budgets], num_worker=num_worker) + data = mix.optimizeMix( + carbon_quota=carbonProdLimit, + time_interval=time_interval, + optimizer=opt_index, + step=sequence, + penalisation=penalisation, + time_index=time_index, + carbon_cost=carbon_cost, + plot="save", + ) + ind_per_opt.update({opt_name: data}) + + for indicator in indicator_list: + new_ind_per_opt = {} + for opt_name, values in ind_per_opt.items(): + ind_per_budget = [] + for budget_value in values: + ind_per_budget.append(float(budget_value[indicator])) + new_ind_per_opt.update({opt_name: ind_per_budget}) + y_tmp.update({indicator: new_ind_per_opt}) + + # plotting + self.plot_evaluation( + X=np.array(budget), + Y=y_tmp, + label_y=indicator_list, + label=optimizer_list, + max_budgets=max_budgets, + average_wide=average_wide, + ) + + # return X, Y, opt_list, max_budgets + return [np.array(budget), y_tmp, optimizer_list, max_budgets] diff --git a/mixsimulator/Experiments/Scenario_type.py b/mixsimulator/Experiments/Scenario_type.py new file mode 100644 index 0000000..a36705b --- /dev/null +++ b/mixsimulator/Experiments/Scenario_type.py @@ -0,0 +1,324 @@ +#!/usr/bin/python3 +import copy +import os +import random +import sys +import threading +import time +import warnings +from datetime import datetime +from math import ceil +from typing import Dict, List + +import matplotlib.pyplot as plt # type: ignore +import numpy as np # type: ignore + +import mixsimulator.nevergradBased.Optimizer as opt +from mixsimulator import ElectricityMix +from mixsimulator.agents.Moderator import StoppableThread + +# from mixsimulator.Evaluation import EvaluationBudget +from mixsimulator.demand.classic.Demand import Demand + +""" + (0) Functions : Check the thread running in background, plot cost results and generate scenarios +""" + + +def generate_random_scenario(centrals: List, time_index: int) -> Dict: + scenario: Dict = {} + tmp: Dict = {} + for central in centrals: + tmp.update({"down": [], "up": []}) + default_proba = random.uniform(0, 0.07) + + for i in range(time_index): + tmp[np.random.choice(["up", "down"], p=[1 - default_proba, default_proba])].append(i) + + up = [] + down = [] + for i in range(time_index): + if i in tmp["down"] and (i - 1 not in tmp["down"] or i == 0): + down.append(i) + if i - 1 in tmp["down"] and i in tmp["up"]: + up.append(i) + tmp["up"] = up + tmp["down"] = down + + scenario.update({central: tmp.copy()}) + + print("scenario: ", scenario) + event_stack: Dict = {} + for i in range(time_index): + for central in scenario.keys(): + if i in scenario[central]["down"]: + try: + event_stack[i].append((central._notify_is_down, central, "down")) + except: + event_stack.update({i: [(central._notify_is_down, central, "down")]}) + elif i in scenario[central]["up"]: + try: + event_stack[i].append((central._notify_is_up, central, "up")) + except: + event_stack.update({i: [(central._notify_is_up, central, "up")]}) + + # print(numpy.arange(0, 2)) + # print("scenario: ", event_stack) + return event_stack + + +def check_thread_running(): + list_ = [] + while True: + tmp = threading.enumerate().copy() + if tmp != list_: + list_ = tmp + for thread in list_: + if thread.is_alive(): + print("THREAD: " + thread.name) + time.sleep(10) + + +def moving_average(x, w): + return np.convolve(x, np.ones(w), "valid") / w + + +def plot_loss(optimum, mode: str = "default", average_wide: int = 0, step: int = 1): + # set the moving average wide + if average_wide == 0: + average_wide = ceil(len(optimum[-1]["coef"]) / 4) + + if mode == "default" or mode == "save" or mode == "best": + # set Y + Y_: Dict[str, List[float]] = {} + Y_["cost"] = [] + Y_["demand_gap"] = [] + for t_i in range(0, len(optimum)): + Y_["cost"].append(optimum[t_i]["loss"]) + Y_["demand_gap"].append(optimum[t_i]["unsatisfied demand"]) + + fig, axs = plt.subplots(1, 2, figsize=(10, 10)) + + # data integration + X = [(i + 1) * step for i in range(0, len(optimum))] + # smooth_value = moving_average(Y_["cost"],average_wide) + for n_axs in range(0, 2): + if n_axs == 0: + # axs[n_axs].plot(X[(average_wide - 1):], smooth_value, '.-' ,alpha=0.5, lw=2, label="cost") + axs[n_axs].plot(X, Y_["cost"], ".-", alpha=0.5, lw=2, label="cost") + else: + axs[n_axs].plot(X, Y_["demand_gap"], "-*", alpha=0.5, lw=2, label="demand gap") + + # Add execution_time and cost information + try: + info = ( + "Final cost: " + + "{:.3f}".format(optimum[-1]["loss"]) + + " ; run_time: " + + "{:.3f}".format(optimum[-1]["elapsed_time"]) + ) + except: + info = "production_cost: " + "{:.3f}".format(optimum[-1]["loss"]) + info += ";Final demand gap: " + "{:.3f}".format(optimum[-1]["unsatisfied demand"]) + plt.annotate( + info, + xy=(0.5, 0), + xytext=(0, 10), + xycoords=("axes fraction", "figure fraction"), + textcoords="offset points", + size=10, + ha="center", + va="bottom", + ) + + # plots parametrizations + for n_axs in range(0, 2): + axs[n_axs].grid() + axs[n_axs].yaxis.set_tick_params(length=0) + axs[n_axs].xaxis.set_tick_params(length=0) + axs[n_axs].set_xlabel("budgets") + # axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) + if n_axs == 0: + axs[n_axs].set_ylabel("cost") + elif n_axs == 1: + axs[n_axs].set_ylabel("Demand gap") + axs[n_axs].legend() + + fig.tight_layout() + try: + path = "cost_results_" + datetime.now().strftime("%d_%m_%Y") + os.makedirs(path) + name = path + "/" + "opt_" + str(datetime.now().strftime("%H%M%S")) + ".png" + fig.savefig(name) + if mode == "default": + plt.show() + + except OSError: + warnings.warn("Can't create the directory " + path + " or already exists") + try: + name = path + "/" + "opt_" + str(datetime.now().strftime("%H%M%S")) + ".png" + fig.savefig(name) + if mode == "default": + plt.show() + except FileNotFoundError: + name = "opt_" + str(datetime.now().strftime("%H%M%S")) + ".png" + fig.savefig(name) + if mode == "default": + plt.show() + + elif mode == "None": + pass + else: + warnings.warn("Choose an available option : default, save or None") + # plt.show() + + +""" +(1) INITIALIZATION PARAMS: + - Configure nevergrad optimizers : argument[1] is opt_name, arg[2] is budget, arg[3] num_worker + - time index as duration (arg[4]), + - step of evaluation as step_budget (arg[5]), + - time_interval as t_i (arg[6]). +""" +opt_name = "OnePlusOne" +budget = 100 +num_worker = 1 +duration = 24 +step_budget = round(budget / 10) +t_i = 1 # time_interval +log_filename = "log.txt" +numb_run = 1 + +for i in range(1, len(sys.argv)): + if i == 1: + opt_name = str(sys.argv[i]) + elif i == 2: + budget = int(sys.argv[i]) + elif i == 3: + num_worker = int(sys.argv[i]) + elif i == 4: + duration = int(sys.argv[i]) + elif i == 5: + step_budget = int(sys.argv[i]) + elif i == 6: + t_i = int(sys.argv[i]) + elif i == 7: + numb_run = int(sys.argv[i]) + elif i == 8: + log_filename = str(sys.argv[i]) + ".pickle" + +for run_ in range(numb_run): + ### optimizer vim + thread_checker = StoppableThread(target=check_thread_running, name="thread_checker") + thread_checker.start() + opt_ = opt.Optimizer(opt=[opt_name], budget=[budget], num_worker=num_worker) + """ + (2) Init MixSimulator instance : + Case one [Default] : "classic" method (see test_classic.py for more use case) + Case two : "MAS" or Multi Agent System method + + Default parameters : + ------------------------ + method : string = "classic", --> method explain above + carbon_cost : float = 0 --> cost of the CO2 production + penalisation_cost: float = 1e7 --> penalisation cost when production is more or less than the demand #NEED VERIFICATION + """ + + classic_mix = ElectricityMix().mix(method="classic", carbon_cost=0, penalisation_cost=100) + mas_mix = ElectricityMix().mix(method="MAS", carbon_cost=0, penalisation_cost=100) + + centrals = mas_mix.get_moderator().get_observable() + scenario = generate_random_scenario(centrals, duration) + + """ + (3) ONE SHOT optimization by calling the classic approach + + """ + classic_mix.set_data_to("Toamasina") + demand = Demand() + demand.set_data_to("Toamasina", delimiter=",") + classic_mix.set_demand(demand) + + start_time = time.time() + classic_result = classic_mix.optimizeMix( + 1e10, optimizer=opt_, step=step_budget, penalisation=100, carbon_cost=0, time_index=duration, plot="save" + ).copy() + ### + ### MODIFY RESULTS BASED ON EVENTS + backup_results = copy.deepcopy(classic_result) + for t in scenario.keys(): + for event in scenario[t]: + for position, central in enumerate(classic_mix.get_centrals()): + if central.get_id() == event[1].get_name(): + if event[2] == "down": + for step_, value in enumerate(classic_result): + for ti in range(t, duration): + classic_result[step_]["coef"][ti][position] = 0.0 + elif event[2] == "up": + for step_, value in enumerate(classic_result): + for ti in range(t, duration): + classic_result[step_]["coef"][ti][position] = backup_results[step_]["coef"][ti][ + position + ] + + for step_, value in enumerate(classic_result): + classic_result[step_]["cost"] = classic_mix.loss_function( + classic_result[step_]["coef"], time_interval=t_i, no_arrange=True + ) + classic_runtime = time.time() - start_time + + ### + ### print(classic_result) + ### PLEASE Check this specific plot in the folder "cost_result_....." + plot_loss(classic_result, step=step_budget, mode="save") + + """ + (4) Simulating the mas platform (Manually) + 1 - First, set params by using set_params method + 2 - Launch the run_optimization method to initiate the simulation + 3 - Add events + 4 - Get the result after all threads done + """ + start_runtime = time.time() + mas_mix.get_moderator().set_params( + 1e10, optimizer=opt_, step=step_budget, penalisation=100, carbon_cost=0, time_index=duration, plot="None" + ) + mas_mix.get_moderator().run_optimization() + + for t in scenario.keys(): + for event in scenario[t]: + event[0](t) + + ### Waiting method + while True: + if len(threading.enumerate()) == 2: + thread_checker.stop() + break + print("SIMULATION DONE") + mas_runtime = time.time() - start_time + + mas_mix.get_moderator().plotResults(mas_mix.get_moderator().get_results(), mode="save") + + ### PLEASE Check this specific plot in the folder "cost_result_....." + plot_loss(mas_mix.get_moderator().get_results(), step=step_budget, mode="save") + + log: Dict = {"run_nb": run_} + log.update( + { + "Perf": ( + (classic_result[-1]["loss"] - mas_mix.get_moderator().get_results()[-1]["loss"]) + / classic_result[-1]["loss"] + ) + * 100, + "events": scenario, + "opt_params": classic_mix.get_params(), + "mas_results": mas_mix.get_moderator().get_results(), + "classic_results": classic_result, + "classic_runtime": classic_runtime, + "mas_runtime": mas_runtime, + } + ) + + file_ = open(log_filename, "a") + print(log, file=file_) + file_.close diff --git a/mixsimulator/LICENSE b/mixsimulator/LICENSE new file mode 100644 index 0000000..bebfcd4 --- /dev/null +++ b/mixsimulator/LICENSE @@ -0,0 +1,21 @@ +MIT License + +Copyright (c) 2020 RASOANAIVO Andry, ANDRIAMALALA Rahamefy Solofohanitra, ANDRIAMIZAKASON Toky Axel + +Permission is hereby granted, free of charge, to any person obtaining a copy +of this software and associated documentation files (the "Software"), to deal +in the Software without restriction, including without limitation the rights +to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +copies of the Software, and to permit persons to whom the Software is +furnished to do so, subject to the following conditions: + +The above copyright notice and this permission notice shall be included in all +copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE +SOFTWARE. diff --git a/mixsimulator/MANIFEST.in b/mixsimulator/MANIFEST.in new file mode 100644 index 0000000..335f07a --- /dev/null +++ b/mixsimulator/MANIFEST.in @@ -0,0 +1,5 @@ +include data/RIToamasina/dataset_RI_Toamasina.csv +include data/RIToamasina/dataset_RI_Toamasina_v2.csv +include data/RIToamasina/dataset_RI_Toamasina_variation_template.csv +include data/RIToamasina/DIR-TOAMASINA_concat.csv +include LICENSE diff --git a/mixsimulator/Mas_platform.py b/mixsimulator/Mas_platform.py new file mode 100644 index 0000000..f2863d4 --- /dev/null +++ b/mixsimulator/Mas_platform.py @@ -0,0 +1,150 @@ +import csv +import pkgutil +import warnings +import tqdm +from typing import Any, List + +import pandas as pd # type: ignore + +from .agents.Moderator import Moderator +from .power_plants.mas.Hydropowerplant import Hydropowerplant +from .power_plants.mas.Thermalpowerplant import Thermalpowerplant + + +class Mas_platform: + """ + Platform of a multi agent system : + - Initialisation of the agents from original dataset + - Get and initiate the optimizer + - Data visualization + """ + + def __init__( + self, + demand=None, + carbon_cost: float = 0, + penalisation_cost: float = 1e7, + bind=None, + raw_data=None, + variation_data=None, + delimiter: str = ";", + variation_delimiter: str = ";", + ) -> None: + self.__moderator = Moderator(carbon_cost, penalisation_cost) + if demand is not None: + self.__moderator.set_demand(demand) + if (bind is not None) and (raw_data is not None): + self.set_data_csv(bind=bind, raw_data=raw_data, delimiter=delimiter) + else: # Load default data + self.set_data_to("Toamasina") + if variation_data is not None: + self.set_variation_csv(bind=variation_data, delimiter=variation_delimiter) + + def get_moderator(self) -> Moderator: + return self.__moderator + + def set_moderator(self, new_moderator) -> None: + self.__moderator = new_moderator + + def set_variation_csv(self, bind=None, delimiter: str = ";") -> None: + if self.__moderator.get_observable() == []: + warnings.warn("Please load a original dataset by using MixSimulator.set_data_csv(...)") + raise + else: + try: + data = pd.DataFrame(pd.read_csv(bind, delimiter=delimiter)) + except FileNotFoundError as e: + print("Error occured on pandas.read_csv : ", e) + print("Please check your file") + raise + except Exception as e: + print("Error occured on pandas.read_csv : ", e) + raise + + for powerplant_index in range(len(self.__moderator.get_observable())): + if self.__moderator.get_observable()[powerplant_index].is_tuneable(): + for i in range(0, data.shape[0]): + if self.__moderator.get_observable()[powerplant_index].get_id() == data["centrals"][i]: + tmp_list: List[Any] = [] + upper: List[Any] = str(data["upper"][i]).split(":") + upper = [float(numeric_string) for numeric_string in upper] + lower: List[Any] = str(data["lower"][i]).split(":") + lower = [float(numeric_string) for numeric_string in lower] + discrete: List[Any] = str(data["discrete"][i]).split(":") + discrete = [float(numeric_string) for numeric_string in discrete] + self.__moderator.get_observable()[powerplant_index].set_variation_params( + lower=lower, upper=upper, choices=discrete + ) + + def set_data_csv(self, bind=None, raw_data=None, delimiter: str = ";"): + if raw_data is not None: + data = pd.DataFrame(raw_data) + # set columns & index + header = data.iloc[0] + data = data[1:] + data.columns = header + data = data.reset_index(drop=True) + for column in data.columns.tolist(): + try: + # convert numeric values + data[column] = pd.to_numeric(data[column]) + except: + pass + else: + try: + data = pd.DataFrame(pd.read_csv(bind, delimiter=delimiter)) + except FileNotFoundError as e: + print("Error occured on pandas.read_csv : ", e) + print("Please check your file") + raise + except Exception as e: + print("Error occured on pandas.read_csv : ", e) + raise + + # self.__reset_centrals() + # powerplant = pc.PowerCentral() + try: + powerplant: Any = ... + self.__moderator.get_demand().set_mean_demand(data["Demand"][0]) + self.__moderator.set_constant_lost(data["lost"][0]) + print("Loading powerplants dataset...") + for i in tqdm.tqdm(range(0, data.shape[0])): + isHydro = data["hydro"][i] + if isHydro.lower() == "true": + powerplant = Hydropowerplant( + data["height"][i], data["flow"][i], data["capacity"][i], data["stock_available"][i], 0.1, 0.8 + ) + else: + powerplant = Thermalpowerplant() + powerplant.set_tuneable(data["tuneable"][i]) + powerplant.set_name(str(data["centrals"][i])) + powerplant.set_fuel_consumption(data["fuel_consumption"][i]) + powerplant.set_availability(data["availability"][i]) + powerplant.set_fuel_cost(data["fuel_cost"][i]) + powerplant.set_initial_value(data["init_value"][i]) + powerplant.set_lifetime(data["lifetime"][i]) + powerplant.set_carbon_prod(data["carbon_production"][i]) + powerplant.set_raw_power(data["raw_power"][i]) + powerplant.set_nb_employees(data["nb_employees"][i]) + powerplant.set_mean_employees_salary(data["mean_salary"][i]) + powerplant.set_max_var(data["max_var"][i]) + ### Register in each agent the observer (moderator) and that allows moderator getting each agent + powerplant.register_observer([self.__moderator]) + except KeyError: + print( + "Columns must be in: tuneable, centrals, fuel_consumption, availability, fuel_cost, init_value, lifetime, carbon_cost, raw_power, nb_employees, mean_salary, demand, lost, height, flow, capacity, stock_available" + ) + raise + + def set_data_to(self, dataset, delimiter: str = ";"): + data: Any = ... + if dataset == "Toamasina": + # by defaut we keep it "Toamasina" + data = pkgutil.get_data("mixsimulator", "/data/RIToamasina/dataset_RI_Toamasina_v2.csv") + data_decoded = csv.reader(data.decode("utf-8").splitlines(), delimiter=delimiter) + self.set_data_csv(raw_data=data_decoded) + else: + # by defaut we keep it "Toamasina" + data = pkgutil.get_data("mixsimulator", "/data/RIToamasina/dataset_RI_Toamasina_v2.csv") + data_decoded = csv.reader(data.decode("utf-8").splitlines(), delimiter=delimiter) + self.set_data_csv(raw_data=data_decoded) diff --git a/mixsimulator/MixSimulator.py b/mixsimulator/MixSimulator.py new file mode 100644 index 0000000..c5fe778 --- /dev/null +++ b/mixsimulator/MixSimulator.py @@ -0,0 +1,503 @@ +import csv +import os +import pkgutil +import warnings +import tqdm +from datetime import datetime +from math import ceil + +# import time +from typing import Any, Dict, List + +import matplotlib.pyplot as plt # type: ignore +import nevergrad as ng +import numpy as np # type: ignore +import pandas as pd # type: ignore + +from .demand.classic.Demand import Demand +from .nevergradBased.Optimizer import Optimizer +from .power_plants.classic import HydroCentral as hc +from .power_plants.classic import PowerCentral as pc + + +class MixSimulator: + """ + The simulator Base + """ + + def __init__(self, carbon_cost: float = 0, penalisation_cost: float = 1000000000000) -> None: + self.__centrals: List[Any] = [] + self.__reset_centrals() + self.__demand = Demand(20, 0.2, 0.2) + self.__lost = 0.0 + self.__penalisation_cost = penalisation_cost + self.__optimizer = Optimizer() + self.__carbon_cost = carbon_cost + self.__carbon_quota = 800139.0 # (g/MWh) + + # for reuse and get_params() + self.time_index = 24 * 7 + self.step = 1 + self.time_interval = 1.0 + self.plot = "default" + self.average_wide = 0 + + def __reset_centrals(self): + self.__centrals = [] + + ## DATA ## + def set_variation_csv(self, bind=None, delimiter: str = ";") -> None: + if self.__centrals == []: + warnings.warn("Please load a original dataset by using MixSimulator.set_data_csv(...)") + raise + else: + try: + data = pd.DataFrame(pd.read_csv(bind, delimiter=delimiter)) + except FileNotFoundError as e: + print("Error occured on pandas.read_csv : ", e) + print("Please check your file") + raise + except Exception as e: + print("Error occured on pandas.read_csv : ", e) + raise + + for central_index in range(len(self.__centrals)): + if self.__centrals[central_index].is_tuneable(): + for i in range(0, data.shape[0]): + if self.__centrals[central_index].get_id() == data["centrals"][i]: + tmp_list: List[Any] = [] + upper: List[Any] = str(data["upper"][i]).split(":") + upper = [float(numeric_string) for numeric_string in upper] + lower: List[Any] = str(data["lower"][i]).split(":") + lower = [float(numeric_string) for numeric_string in lower] + choices = None + if str(data["discrete"][i]) != "None": + discrete: List[Any] = str(data["discrete"][i]).split(":") + discrete = [float(numeric_string) for numeric_string in discrete] + choices = discrete + self.__centrals[central_index].set_variation_params( + lower=lower, upper=upper, choices=choices + ) + + def set_data_csv(self, bind=None, raw_data=None, delimiter: str = ";"): + if raw_data is not None: + data = pd.DataFrame(raw_data) + # set columns & index + header = data.iloc[0] + data = data[1:] + data.columns = header + data = data.reset_index(drop=True) + for column in data.columns.tolist(): + try: + # convert numeric values + data[column] = pd.to_numeric(data[column]) + except: + pass + else: + try: + data = pd.DataFrame(pd.read_csv(bind, delimiter=delimiter)) + except FileNotFoundError as e: + print("Error occured on pandas.read_csv : ", e) + print("Please check your file") + raise + except Exception as e: + print("Error occured on pandas.read_csv : ", e) + raise + + self.__reset_centrals() + centrale = pc.PowerCentral() + try: + print("Loading powerplants dataset...") + for i in tqdm.tqdm(range(0, data.shape[0])): + isHydro = data["hydro"][i] + if isHydro == True: + centrale = hc.HydroCentral( + data["height"][i], data["flow"][i], data["capacity"][i], data["stock_available"][i], 0.1, 0.8 + ) + else: + centrale = pc.PowerCentral() + centrale.set_tuneable(data["tuneable"][i]) + centrale.set_id(str(data["centrals"][i])) + centrale.set_fuel_consumption(data["fuel_consumption"][i]) + centrale.set_availability(data["availability"][i]) + centrale.set_fuel_cost(data["fuel_cost"][i]) + centrale.set_initial_value(data["init_value"][i]) + centrale.set_lifetime(data["lifetime"][i]) + centrale.set_carbon_prod(data["carbon_production"][i]) + centrale.set_raw_power(data["raw_power"][i]) + centrale.set_nb_employees(data["nb_employees"][i]) + centrale.set_mean_employees_salary(data["mean_salary"][i]) + centrale.set_max_var(data["max_var"][i]) + self.__centrals.append(centrale) + self.__demand.set_mean_demand(data["Demand"][0]) + self.__lost = data["lost"][0] + except KeyError: + print( + "Columns must be in: tuneable, centrals, fuel_consumption, availability, fuel_cost, init_value, lifetime, carbon_cost, raw_power, nb_employees, mean_salary, demand, lost, height, flow, capacity, stock_available" + ) + raise + + def set_data_to(self, dataset, delimiter: str = ";"): + data: Any = ... + if dataset == "Toamasina": + # by defaut we keep it "Toamasina" + data = pkgutil.get_data("mixsimulator", "/data/RIToamasina/dataset_RI_Toamasina_v2.csv") + data_decoded = csv.reader(data.decode("utf-8").splitlines(), delimiter=delimiter) + self.set_data_csv(raw_data=data_decoded) + else: + # by defaut we keep it "Toamasina" + data = pkgutil.get_data("mixsimulator", "/data/RIToamasina/dataset_RI_Toamasina_v2.csv") + data_decoded = csv.reader(data.decode("utf-8").splitlines(), delimiter=delimiter) + self.set_data_csv(raw_data=data_decoded) + + def set_demand(self, demand: Demand): + self.__demand = demand + + def get_demand(self): + return self.__demand + + def get_centrals(self): + return self.__centrals + + def set_lost(self, lost: float): + self.__lost = lost + + def set_optimizer(self, optimizer: Optimizer): + self.__optimizer = optimizer + + def get_optimizer(self) -> Optimizer: + return self.__optimizer + + def set_carbon_quota(self, cb_quota: float): + self.__carbon_quota = cb_quota + + def set_penalisation_cost(self, k): + self.__penalisation_cost = k + + def get_penalisation_cost(self): + return self.__penalisation_cost + + def set_carbon_cost(self, carbon_cost): + self.__carbon_cost = carbon_cost + + def get_carbon_cost(self): + return self.__carbon_cost + + ## EVALUATION FONCTION ## + def get_production_cost_at_t(self, usage_coef, time_index, time_interval): + production_cost = 0 + for centrale_index in range(0, len(self.__centrals)): + central = self.__centrals[centrale_index] + production_cost += ( + ( + central.get_fuel_cost() + * central.get_fuel_consumption() + * central.get_raw_power() + * usage_coef[centrale_index] + ) + + (central.get_employees_salary() + central.get_amortized_cost(time_index)) + ) * time_interval + return production_cost + + def get_production_at_t(self, usage_coef, time_interval): + production = 0 + for centrale_index in range(0, len(self.__centrals)): + central = self.__centrals[centrale_index] + production += (central.get_raw_power() * usage_coef[centrale_index]) * time_interval + return production + + def get_unsatisfied_demand_at_t(self, usage_coef, time_index, time_interval): + # return ( self.__demand.get_demand_approxima(time_index, time_interval) - self.get_production_at_t(usage_coef, time_interval)) + return self.__demand.get_demand_monthly(time_index, time_interval) - self.get_production_at_t( + usage_coef, time_interval + ) + + def get_carbon_production_at_t(self, usage_coef, time_interval): + carbon_prod = 0 + for centrale_index in range(0, len(self.__centrals)): + central = self.__centrals[centrale_index] + carbon_prod += ( + central.get_carbon_production() * usage_coef[centrale_index] * central.get_raw_power() + ) * time_interval + return carbon_prod + + def get_carbon_over_production(self, usage_coef, time_interval): + emited_carbon = 0 # (g) + total_production = 1 + for t in range(0, len(usage_coef)): + emited_carbon += self.get_carbon_production_at_t(usage_coef[t], time_interval) + total_production += self.get_production_at_t(usage_coef[t], time_interval) + carbon_production = emited_carbon / total_production + return max(0, carbon_production - self.__carbon_quota) # (g/MWh) + + def get_weighted_coef(self, usage_coef, time_interval): + for central_index in range(0, len(usage_coef[0])): + self.__centrals[central_index].reset_central() + weighted_coef = usage_coef.copy() + for t in range(0, len(weighted_coef)): + for central_index in range(0, len(weighted_coef[t])): + min_av = self.__centrals[central_index].get_min_availability(t) + max_av = self.__centrals[central_index].get_max_availability(t) + if max_av < min_av: + min_av = 0 # a verifier + weighted_coef[t][central_index] = min_av + weighted_coef[t][central_index] * (max_av - min_av) + self.__centrals[central_index].back_propagate(weighted_coef[t][central_index], t, time_interval) + return weighted_coef + + def loss_function(self, usage_coef, time_interval: int = 1, no_arrange=False, init: int = 0) -> float: + if no_arrange is False: + usage_coef = self.arrange_coef_as_array_of_array(usage_coef) + weighted_coef = self.get_weighted_coef(usage_coef, time_interval=time_interval) + else: + weighted_coef = usage_coef + loss = 0 + for t in range(0, len(weighted_coef)): + loss += self.get_production_cost_at_t(weighted_coef[t], t, time_interval) + ( + self.get_penalisation_cost() + * np.abs(self.get_unsatisfied_demand_at_t(weighted_coef[t], t, time_interval)) + ) + loss += self.get_carbon_cost() * (self.get_carbon_over_production(weighted_coef, time_interval)) + return loss + + def arrange_coef_as_array_of_array(self, raw_usage_coef): + ordered_coef = [] + cur_time_coef = [] + for coef_index in range(len(raw_usage_coef)): + cur_time_coef.append(raw_usage_coef[coef_index]) + ## indice de la premiere centrale a t+1 + if (coef_index + 1) % len(self.__centrals) == 0: + ordered_coef.append(cur_time_coef) + cur_time_coef = [] + return ordered_coef + + ## OPTiMiZATION ## + + def __opt_params(self, time_index): + self.__optimizer.set_parametrization(self.get_opt_params(time_index)) + + def get_opt_params(self, time_index): + variable_parametrization = [] + for _ in range(time_index): + for central_index in range(len(self.__centrals)): + if not self.__centrals[central_index].is_tuneable(): + variable_parametrization += [ng.p.Choice([0.0, 1.0])] + else: + # check the params by + # print(self.__centrals[central_index].get_variation_params()) + variable_parametrization += [self.__centrals[central_index].get_variation_params()] + return ng.p.Tuple(*variable_parametrization) + + def optimizeMix( + self, + carbon_quota: float = None, + demand: Demand = None, + lost: float = None, + optimizer: Optimizer = None, + step: int = 1, + time_index: int = 24 * 7, + time_interval: float = 1, + penalisation: float = None, + carbon_cost: float = None, + plot: str = "default", + average_wide: int = 0, + ): + + # init params + self.time_index = time_index + self.step = step + self.time_interval = time_interval + self.plot = plot + self.average_wide = average_wide + if demand is not None: + self.set_demand(demand) + if lost is not None: + self.set_lost(lost) + if penalisation is not None: + self.set_penalisation_cost(penalisation) + if carbon_cost is not None: + self.set_carbon_cost(carbon_cost) + if carbon_quota is not None: + self.set_carbon_quota(carbon_quota) + if optimizer is not None: + self.set_optimizer(optimizer) + + # tune optimizer parametrization + self.__opt_params(self.time_index) + + # init constraints + constraints = {} + constraints.update({"time_interval": self.time_interval}) + + # let's optimize + results = self.__optimizer.optimize( + mix=self, func_to_optimize=self.loss_function, constraints=constraints, step=self.step + ) + + results = self.__reshape_results(results, self.time_interval) + + self.plotResults(results, mode=self.plot, time_interval=self.time_interval, average_wide=self.average_wide) + + return results + + def __reshape_results(self, results, time_interval): + for tmp in results: + usage_coef = self.arrange_coef_as_array_of_array(tmp["coef"]) + tmp.update({"coef": self.get_weighted_coef(usage_coef, time_interval)}) + + # for central_index in range(0, len(self.__centrals)): + # try: + # print(self.__centrals[central_index].get_stock_evolution()) + # # Not a hydro power plant, so the methode does not exist + # except: + # pass + return results + + def get_params(self) -> dict: + return { + "centrals": self.__centrals, + "optimizer": self.get_optimizer(), + "penalisation_cost": self.get_penalisation_cost(), + "carbon_cost": self.get_carbon_cost(), + "demand": self.__demand, + "lost": self.__lost, + "carbon_quota": self.__carbon_quota, + "step": self.step, + "time_interval": self.time_interval, + "time_index": self.time_index, + "plot": self.plot, + "moving average_wide": self.average_wide, + } + + ## PLOT ## + def moving_average(self, x, w): + return np.convolve(x, np.ones(w), "valid") / w + + def plotResults(self, optimum: dict = {}, mode: str = "default", time_interval: float = 1, average_wide: int = 0): + # set the moving average wide + if average_wide == 0: + average_wide = ceil(len(optimum[-1]["coef"]) / 4) + + if mode == "default" or mode == "save" or mode == "best": + # set Y of axe 1 + Y: Dict[str, List[float]] = {} + label_y: List[str] = [] + sum_prod = [] + for c in self.__centrals: + label_y.append(c.get_id()) + Y.update({c.get_id(): []}) + if mode != "best": + for array in optimum[-1]["coef"]: + sum_t = 0.0 + for enum, value in enumerate(array): + Y[label_y[enum]].append(value) + sum_t = sum_t + (value * self.__centrals[enum].get_raw_power()) + sum_prod.append(sum_t) + else: + ##TODO get best loss + pass + + # set Y of axe 2 + Y_: Dict[str, List[float]] = {} + Y_["demand"] = [] + Y_["production"] = sum_prod + for t_i in range(0, len(optimum[-1]["coef"])): + Y_["demand"].append(self.__demand.get_demand_monthly(t_i, time_interval)) + + fig, axs = plt.subplots(1, 2, figsize=(12, 6)) + + # data integration + X = [i for i in range(0, len(optimum[-1]["coef"]))] + for n_axs in range(0, 2): + if n_axs == 0: + for central, value in Y.items(): + smooth_value = self.moving_average(value, average_wide) + axs[n_axs].plot(X[(average_wide - 1) :], smooth_value, ".-", alpha=0.5, lw=2, label=central) + elif n_axs == 1: + axs[n_axs].plot(X, Y_["demand"], label="demand") + smooth_value = self.moving_average(Y_["production"], average_wide) + axs[n_axs].plot(X[(average_wide - 1) :], smooth_value, ".-", alpha=0.5, lw=2, label="production") + + # Add execution_time and loss information + try: + info = ( + "production_cost: " + + "{:.3f}".format(optimum[-1]["loss"]) + + " ; execution_time: " + + "{:.3f}".format(optimum[-1]["elapsed_time"]) + ) + except: + info = "production_cost: " + "{:.3f}".format(optimum[-1]["loss"]) + info += "; demand gap: " + "{:.3f}".format(optimum[-1]["unsatisfied demand"]) + plt.annotate( + info, + xy=(0.5, 0), + xytext=(0, 10), + xycoords=("axes fraction", "figure fraction"), + textcoords="offset points", + size=10, + ha="center", + va="bottom", + ) + + # plots parametrizations + for n_axs in range(0, 2): + axs[n_axs].grid() + axs[n_axs].yaxis.set_tick_params(length=0) + axs[n_axs].xaxis.set_tick_params(length=0) + axs[n_axs].set_xlabel("hours") + # axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) + if n_axs == 0: + axs[n_axs].set_ylabel("coef (%)") + elif n_axs == 1: + axs[n_axs].set_ylabel("MWh") + axs[n_axs].legend() + + fig.tight_layout() + try: + path = "results_" + datetime.now().strftime("%d_%m_%Y") + os.makedirs(path) + name = ( + path + + "/" + + "opt_" + + str(self.get_optimizer().get_optimizers()) + + "_" + + datetime.now().strftime("%H%M%S") + + ".png" + ) + fig.savefig(name) + if mode == "default": + plt.show() + + except OSError: + warnings.warn("Can't create the directory " + path + " or already exists") + try: + name = ( + path + + "/" + + "opt_" + + str(self.get_optimizer().get_optimizers()) + + "_" + + datetime.now().strftime("%H%M%S") + + ".png" + ) + fig.savefig(name) + if mode == "default": + plt.show() + except FileNotFoundError: + name = ( + "opt_" + + str(self.get_optimizer().get_optimizers()) + + "_" + + datetime.now().strftime("%H%M%S") + + ".png" + ) + fig.savefig(name) + if mode == "default": + plt.show() + + elif mode == "None": + pass + else: + warnings.warn("Choose an available option : default, save or None") + # plt.show() \ No newline at end of file diff --git a/mixsimulator/README.md b/mixsimulator/README.md new file mode 100644 index 0000000..396fab3 --- /dev/null +++ b/mixsimulator/README.md @@ -0,0 +1,77 @@ +# MixSimulator +MixSimulator is an application with an optimization model for calculating and simulating the least cost of an energy mix under certain constraints. The optimizers used are based on the Nevergrad Python 3.6+ library. + +The primary objective of the simulator is to study the relevance of an energy mix connected to each Inter-connected Grid through the coefficient of usage of each unit in the production cost. + +## Version 0.4 +The current version is a multi-agent system (MAS) approach but keeps the previous classic optimization approach available. Check `test_mas.py` for more details. (Available events are : powerplant shutdown and powerplant power_on). + +An experiment on evaluating both method is available in `Experiments/Scenario_type.py` or [![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/drive/1IzjJlNg3fCp14GICB2GGYYwmSIEK9yQf?usp=sharing) + +#### Note +This project is a work in progress so it can not yet used in production (Many changes are on their way). Feedbacks are welcome! + +### Specifications : +- Generic simulator, compatible with data from Madagascar and those from abroad (but may require data pre-processing beforehand); +- Optimization of duty cycle (or usage coefficient) under constraints ; +- Get the production cost and various performance indicators (CO2 emission, unsatisfied demand); +- Estimate of the costs of a mix or a power plant over the long term ; +- Comparison of the performance indicators on different optimizers. (see `EvaluationBudget` method) + +### Perspectives : +- Add other constraints ; +- Long-term Optimization ; +- Pair with a transmission and distribution power grid simulator (MixSimulator can provide input data). + +Suggestions are welcome! + +## How to install +It can be installed with : +```python +pip install mixsimulator +``` +MixSimulator is written in Python 3.6 and requires the following Python packages : nevergrad, prophet, typing, numpy, pandas and matplotlib. + +## How to run +As MixSimulator is a python package, it can be called and used as we can see in `test_classic.py` and `test_mas.py`. + +Official documentation will accompany the first release version. + +## Datasets + +### Power plants dataset +The dataset "dataset_RI_Toamasina_v2.csv" is for the test and it comes from the Inter-connected energy mix of Toamasina Madagascar (2017) and Some information from the dataset is estimated. + +Dataset features needed: +- `tuneable` (boolean): is the control unit controllable or not? +- `green` (boolean): is it a renewable energy plant? +- `hydro` (boolean): is it a hydro power plant? +- `fuel` (boolean): is it a thermal power plant? +- `centrals` : identification +- `fuel_consumption` (g/MWh): fuel consumption (in the case of a fossil fuel power plant) +- `availability` (%): plant availability +- `fuel_cost` ($/g): price of fuel used +- `init_value` ($): initial investment in setting up the plant +- `lifetime` (years): plant lifetime +- `carbon_production` (g/MWh): emission rate of CO2 from the power plant +- `raw_power` (MW): nominal power of the plant +- `nb_employees`: number of employees at the central level +- `mean_salary` ($): average salary of plant employees +- `demand` (MWh): electricity demand +- `lost` (MWh): electrical loss at another level (ie: transport network) + +Hydro specification : +- `height` (meter): height of the stream ; +- `flow` : flow of the stream ; +- `stock_available` : water reservoir ; +- `capacity` : max water reservoir. + +`nb_employees * mean_salary` **can be used as a variable cost of the plant if you want to directly use other informations as variable cost.** + +### Demand and Variation datas +There is also "DIR-TOAMASINA_concat.csv" about Consumption data (in kwh, more details in `demand/`) and "dataset_RI_Toamasina_variation_template.csv" about limits in variation of power plants load following (WIP). + +**If you have datasets of any region in the world that can be used to evaluate our model, please contact us.** + +## Contact +For questions and feedbacks related to the project, please send an email to r.andry.rasoanaivo@gmail.com or soloforahamefy@gmail.com or tokyandriaxel@gmail.com diff --git a/mixsimulator/__init__.py b/mixsimulator/__init__.py new file mode 100644 index 0000000..8f0f0a5 --- /dev/null +++ b/mixsimulator/__init__.py @@ -0,0 +1,16 @@ +# This source code is licensed under the MIT license found in the +# LICENSE file in the root directory of this source tree. + +from .base import ElectricityMix +from .Mas_platform import Mas_platform +from .MixSimulator import MixSimulator + +# from .Demand import Demand +from .nevergradBased import Optimizer +from .Evaluation import EvaluationBudget + +# __all__ = ["ElectricityMix","Mas_platform","MixSimulator","Demand","Optimizer","EvaluationBudget"] +__all__ = ["ElectricityMix", "Mas_platform", "MixSimulator", "Optimizer", "EvaluationBudget"] + + +__version__ = "0.4" diff --git a/mixsimulator/agents/Agent.py b/mixsimulator/agents/Agent.py new file mode 100644 index 0000000..73706ed --- /dev/null +++ b/mixsimulator/agents/Agent.py @@ -0,0 +1,88 @@ +import json +import pkgutil +import uuid +from threading import Lock, Timer +from typing import Any, Dict, List + +from .Interfaces import Observable + + +class Periodic: + """ + A periodic task running in threading.Timers + """ + + def __init__(self): + self._lock = Lock() + self._timer = None + self._stopped = True + + def set_function(self, function, *args, **kwargs): + self.__function = function + self.args = args + self.kwargs = kwargs + + def set_interval(self, interval): + self.__interval = interval + + def start(self, from_run=False): + self._lock.acquire() + if from_run or self._stopped: + self._stopped = False + self._timer = Timer(self.__interval, self._run) + self._timer.start() + self._lock.release() + + def _run(self): + self.start(from_run=True) + self.__function(*self.args, **self.kwargs) + + def stop(self): + self._lock.acquire() + self._stopped = True + self._timer.cancel() + self._lock.release() + + +class Agent(Observable): + def __init__(self) -> None: + super().__init__() + self._code_files = "params_files/exchange_code.json" + self.__type = "empty" + self.__name = "" + self._scheduled_actions: Dict[Any, Any] = {} + self.set_id(uuid.uuid4()) + + def __repr__(self): + return self.get_name() + "[" + str(self.get_id()) + "] (" + self.__type + ")" + + def _schedule_action(self, actions: Dict) -> None: + ## action is a dict {function, interval} + for action in actions.keys(): + if action not in self._scheduled_actions.keys(): + ## self._schedule_actions is a dict {function, Periodic --the class above--)} + periodic = Periodic() + periodic.set_function(action) + periodic.set_interval(actions[action]) + self._scheduled_actions.update({action: periodic}) + self._scheduled_actions[action].start() + + def register_observer(self, moderators: List, t_from=0) -> None: + super().register_observer(moderators) + data_json: Any = pkgutil.get_data("mixsimulator", self._code_files) + register_signal: Any = json.loads(data_json.decode("utf-8"))["100"] + register_signal["id"] = self.get_id() + register_signal["t_from"] = t_from + self._notify_observers(register_signal) + + def set_type(self, ntype): + self.__type = ntype + + def get_type(self): + return self.__type + + def set_name(self, name): + self.__name = name + + def get_name(self): + return self.__name diff --git a/mixsimulator/agents/Interfaces.py b/mixsimulator/agents/Interfaces.py new file mode 100644 index 0000000..3e3bac3 --- /dev/null +++ b/mixsimulator/agents/Interfaces.py @@ -0,0 +1,30 @@ +from typing import Any, List + + +class Observable: + def __init__(self): + self.__observers = [] + self.__id: Any = None + + def get_id(self) -> str: + return self.__id + + def set_id(self, id) -> None: + self.__id = id + + def get_observers(self) -> List: + return self.__observers + + def register_observer(self, observers: List) -> None: + for observer in observers: + if observer not in self.__observers: + self.__observers.append(observer) + + def _notify_observers(self, signal) -> None: + for obs in self.get_observers(): + obs.update(self, signal) + + +class Observer: + def update(self, observable, *args, **kwargs): + pass diff --git a/mixsimulator/agents/Moderator.py b/mixsimulator/agents/Moderator.py new file mode 100644 index 0000000..52ffa9a --- /dev/null +++ b/mixsimulator/agents/Moderator.py @@ -0,0 +1,630 @@ +import copy +import ctypes +import os +import threading +import warnings +from datetime import datetime +from math import ceil + +# import time +from typing import Any, Dict, List + +import matplotlib.pyplot as plt # type: ignore +import nevergrad as ng +import numpy as np # type: ignore + +from ..demand.mas.Demand import Demand +from ..nevergradBased.Optimizer import Optimizer +from ..power_plants.mas.PowerPlant import PowerPlant +from .Interfaces import Observable, Observer + + +class StoppableThread(threading.Thread): + """ + Thread class with a stop() method. The thread itself has to check + regularly for the stopped() condition. + """ + + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) + + def run(self): + try: + super().run() + finally: + # print('stopped thread: ' + self.name) + pass + + def get_id(self): + if hasattr(self, "_thread_id"): + return self._thread_id + for id, thread in threading._active.items(): + if thread is self: + return id + + def stop(self): + thread_id = self.get_id() + resu = ctypes.pythonapi.PyThreadState_SetAsyncExc(ctypes.c_long(thread_id), ctypes.py_object(SystemExit)) + if resu > 1: + ctypes.pythonapi.PyThreadState_SetAsyncExc(thread_id, 0) + + +class Moderator(Observer): + def __init__(self, carbon_cost, penalisation_cost) -> None: + super().__init__() + self.__log_file = open("optimisation.log", "w") + self.__observable: List[PowerPlant] = [] + self.__powerplants_down: List[PowerPlant] = [] + dm = Demand(demand=20, var_per_day=0.2, var_per_season=0.2) + dm.set_data_to("Toamasina", delimiter=",") + self.__params_state: Any = None + self.__demand = dm + self.__cst_lost = 0.0 + self.__penalisation_cost = penalisation_cost + self.__optimizer = Optimizer() + self.__carbon_cost = carbon_cost + self.__carbon_quota = 800139.0 # (g/MWh) + self.__results: List[Any] = [] + self.__latest_results: List[Any] = [] + + # for reuse and get_params() + self.time_index = 24 * 7 + self.step = 1 + self.time_interval = 1.0 + self.plot = "default" + self.average_wide = 0 + + # self.__cur_optimization: StoppableThread = None + self.__cur_optimization: Any = None + + def __reset_powerplant(self): + self.__observable: List[PowerPlant] = [] + + def get_observable(self) -> List[PowerPlant]: + return self.__observable + + def set_observable(self, observables: List[PowerPlant]) -> None: + self.__observable = observables + + def __add_observable(self, observable: PowerPlant) -> None: + if observable not in self.__observable: + self.__observable.append(observable) + + def __shutdown_observable(self, observable: PowerPlant, id) -> None: + observable.shutdown() + + def __power_on_observable(self, observable: PowerPlant) -> None: + observable.power_on() + + ### COMMUNICATION + def __update_self(self, observable, signal): + print(observable, "sends signal code ", signal["code"]) + if signal["code"] == 100: + self.__add_observable(observable) + elif signal["code"] == 400: + self.__shutdown_observable(observable, signal["id"]) + elif signal["code"] == 200: + self.__power_on_observable(observable) + + def run_optimization(self): + """ + Initial run of the simulation (must be run at first) + """ + _kwargs = { + "carbon_quota": self.__carbon_quota, + "demand": self.get_demand(), + "lost": self.__cst_lost, + "optimizer": self.__optimizer, + "step": self.step, + "time_index": self.time_index, + "time_interval": self.time_interval, + "penalisation": self.__penalisation_cost, + "carbon_cost": self.get_carbon_cost(), + "plot": self.plot, + "average_wide": self.average_wide, + } + self.__cur: int = 0 + self.__cur_optimization = StoppableThread(target=self.optimizeMix, kwargs=_kwargs, name="Initial run") + self.__cur_optimization.start() + + def update(self, observable, signal, *args, **kwargs) -> None: # type: ignore + if self.__params_state is None: + self.__update_self(observable, signal) + else: + ## relaunching optimization + _kwargs = { + "carbon_quota": self.__carbon_quota, + "demand": self.get_demand(), + "lost": self.__cst_lost, + "optimizer": self.__optimizer, + "step": self.step, + "time_index": self.time_index, + "time_interval": self.time_interval, + "penalisation": self.__penalisation_cost, + "carbon_cost": self.get_carbon_cost(), + "plot": self.plot, + "average_wide": self.average_wide, + } + if self.__cur_optimization is None: + pass + elif signal["t_from"] == self.__cur: + _kwargs.update({"init": signal["t_from"], "time_index": _kwargs["time_index"] - signal["t_from"]}) + self.__cur_optimization.stop() + self.__update_self(observable, signal) + self.__cur_optimization = StoppableThread( + target=self.optimizeMix, kwargs=_kwargs, name=str(observable) + str(signal["code"]) + ) + self.__cur = signal["t_from"] + self.__cur_optimization.start() + elif signal["t_from"] > self.__cur: + _kwargs.update({"init": signal["t_from"], "time_index": _kwargs["time_index"] - signal["t_from"]}) + print("wainting for current optimization to finish".upper()) + self.__cur_optimization.join() + self.__update_self(observable, signal) + self.__cur_optimization = StoppableThread( + target=self.optimizeMix, kwargs=_kwargs, name=str(observable) + str(signal["code"]) + ) + self.__cur = signal["t_from"] + self.__cur_optimization.start() + + ### PARAMETRIZATION + def set_demand(self, demand_agent) -> None: + self.__demand = demand_agent + + def get_demand(self): + return self.__demand + + def set_mean_demand(self, mean_demand) -> None: + self.__demand = mean_demand + + def set_constant_lost(self, lost) -> None: + self.__cst_lost = lost + + def set_optimizer(self, optimizer: Optimizer): + self.__optimizer = optimizer + + def get_optimizer(self) -> Optimizer: + return self.__optimizer + + def set_carbon_quota(self, cb_quota: float): + self.__carbon_quota = cb_quota + + def set_penalisation_cost(self, k): + self.__penalisation_cost = k + + def get_penalisation_cost(self): + return self.__penalisation_cost + + def set_carbon_cost(self, carbon_cost): + self.__carbon_cost = carbon_cost + + def get_carbon_cost(self): + return self.__carbon_cost + + def get_results(self): + return self.__results + + ### EVALUATION FUNCTIONS + def get_production_cost_at_t(self, usage_coef, time_index, time_interval): + production_cost = 0 + for powerplant_index in range(0, len(self.__observable)): + if self.__observable[powerplant_index].get_type() != "Demand": + powerplant = self.__observable[powerplant_index] + production_cost += ( + ( + powerplant.get_fuel_cost() + * powerplant.get_fuel_consumption() + * powerplant.get_raw_power() + * usage_coef[powerplant_index] + ) + + (powerplant.get_employees_salary() + powerplant.get_amortized_cost(time_index)) + ) * time_interval + return production_cost + + def get_production_at_t(self, usage_coef, time_interval): + production = 0 + for powerplant_index in range(0, len(self.__observable)): + if self.__observable[powerplant_index].get_type() != "Demand": + powerplant = self.__observable[powerplant_index] + production += (powerplant.get_raw_power() * usage_coef[powerplant_index]) * time_interval + return production + + def get_unsatisfied_demand_at_t(self, usage_coef, time_index, time_interval, init: int = 0): + # return ( self.__demand.get_demand_approxima(time_index, time_interval) - self.get_production_at_t(usage_coef, time_interval)) + return self.__demand.get_demand_monthly(time_index, time_interval, init=init) - self.get_production_at_t( + usage_coef, time_interval + ) + + def get_carbon_production_at_t(self, usage_coef, time_interval): + carbon_prod = 0 + for powerplant_index in range(0, len(self.__observable)): + if self.__observable[powerplant_index].get_type() != "Demand": + powerplant = self.__observable[powerplant_index] + carbon_prod += ( + powerplant.get_carbon_production() * usage_coef[powerplant_index] * powerplant.get_raw_power() + ) * time_interval + return carbon_prod + + def get_carbon_over_production(self, usage_coef, time_interval): + emited_carbon = 0 # (g) + total_production = 1 + for t in range(0, len(usage_coef)): + emited_carbon += self.get_carbon_production_at_t(usage_coef[t], time_interval) + total_production += self.get_production_at_t(usage_coef[t], time_interval) + carbon_production = emited_carbon / total_production + return max(0, carbon_production - self.__carbon_quota) # (g/MWh) + + def get_weighted_coef(self, usage_coef, time_interval, init: int = 0): + if init == 0: + for powerplant in self.__observable: + if powerplant.get_type() not in ["Demand"]: + powerplant.reset_powerplant() + else: + for powerplant in self.__observable: ########################################### + if powerplant.get_type() not in [ + "Demand", + "Hydropowerplant", + ]: ## Keep previous history stock for hydro ## + powerplant.reset_powerplant() ########################################### + + weighted_coef = usage_coef.copy() + for t in range(0, len(weighted_coef)): + for powerplant_index, powerplant in enumerate(self.__observable): + if powerplant.get_type() != "Demand": + min_av = powerplant.get_min_availability() + max_av = powerplant.get_max_availability(t, init=init) + if max_av < min_av: + min_av = 0 # a verifier + weighted_coef[t][powerplant_index] = min_av + weighted_coef[t][powerplant_index] * ( + max_av - min_av + ) + powerplant.back_propagate(weighted_coef[t][powerplant_index], t, time_interval) + return weighted_coef + + def loss_function(self, usage_coef, time_interval: int = 1, no_arrange=False, init: int = 0) -> float: + if no_arrange is False: + usage_coef = self.arrange_coef_as_array_of_array(usage_coef) + weighted_coef = self.get_weighted_coef(usage_coef, time_interval=time_interval, init=init) + else: + weighted_coef = usage_coef + loss = 0 + for t in range(0, len(weighted_coef)): + loss += self.get_production_cost_at_t(weighted_coef[t], t, time_interval) + ( + self.get_penalisation_cost() + * np.abs(self.get_unsatisfied_demand_at_t(weighted_coef[t], t, time_interval, init=init)) + ) + loss += self.get_carbon_cost() * (self.get_carbon_over_production(weighted_coef, time_interval)) + return loss + + def arrange_coef_as_array_of_array(self, raw_usage_coef): + ordered_coef = [] + cur_time_coef = [] + for coef_index in range(len(raw_usage_coef)): + cur_time_coef.append(raw_usage_coef[coef_index]) + ## indice de la premiere powerplant a t+1 + if (coef_index + 1) % len(self.__observable) == 0: + ordered_coef.append(cur_time_coef) + cur_time_coef = [] + return ordered_coef + + ### OPTiMiZATION + def __opt_params(self, time_index): + self.__optimizer.set_parametrization(self.get_opt_params(time_index)) + + def get_opt_params(self, time_index): + variable_parametrization = [] + for _ in range(time_index): + for powerplant_index in range(len(self.__observable)): + if self.__observable[powerplant_index].get_type() != "Demand": + if not self.__observable[powerplant_index].is_tuneable(): + variable_parametrization += [ng.p.Choice([0.0, 1.0])] + else: + # check the params by + variable_parametrization += [self.__observable[powerplant_index].get_variation_params()] + return ng.p.Tuple(*variable_parametrization) + + def __reshape_results(self, results, time_interval): + for tmp in results: + usage_coef = self.arrange_coef_as_array_of_array(tmp["coef"]) + tmp.update({"coef": self.get_weighted_coef(usage_coef, time_interval)}) + + # for powerplant_index in range(0, len(self.__observable)): + # try: + # print(self.__observable[powerplant_index].get_stock_evolution()) + # # Not a hydro power plant, so the methode does not exist + # except: + # pass + return results + + def optimizeMix( + self, + carbon_quota: float = None, + demand: Demand = None, + lost: float = None, + optimizer: Optimizer = None, + step: int = 1, + time_index: int = 24 * 7, + time_interval: float = 1, + penalisation: float = None, + carbon_cost: float = None, + plot: str = "default", + average_wide: int = 0, + init: int = 0, + ): + + # init params + print( + "OPTIMIZATION IS RUNNING:\n time_index = " + + str(time_index) + + " hours\n from init=" + + str(init) + + " to " + + str(self.time_index) + ) + # step is the step of opt budget + if time_index < time_interval: + time_interval = time_index + else: + time_interval = time_interval + if demand is not None: + self.set_demand(demand) + if lost is not None: + self.set_constant_lost(lost) + if penalisation is not None: + self.set_penalisation_cost(penalisation) + if carbon_cost is not None: + self.set_carbon_cost(carbon_cost) + if carbon_quota is not None: + self.set_carbon_quota(carbon_quota) + if optimizer is not None: + self.set_optimizer(optimizer) + + # tune optimizer parametrization + self.__opt_params(time_index) + + # init constraints + constraints = {} + constraints.update({"time_interval": time_interval}) + + # let's optimize + self.__latest_results = self.__optimizer.optimize( + mix=self, func_to_optimize=self.loss_function, constraints=constraints, step=step, init=init + ) + + self.__latest_results = self.__reshape_results(self.__latest_results, time_interval) + + self.plotResults(self.__latest_results, mode="save", time_interval=time_interval, average_wide=average_wide) + + # update results when init is different of 0 + if init == 0: + self.__results = copy.deepcopy(self.__latest_results) + else: + self.__results = self.__update_results(self.__results, self.__latest_results, init) + + print("optimization is done".upper()) + return self.__results + + def get_params(self) -> dict: + return { + "agents": self.__observable, + "optimizer": self.get_optimizer(), + "penalisation_cost": self.get_penalisation_cost(), + "carbon_cost": self.get_carbon_cost(), + "demand": self.__demand, + "lost": self.__cst_lost, + "carbon_quota": self.__carbon_quota, + "step": self.step, + "time_interval": self.time_interval, + "time_index": self.time_index, + "plot": self.plot, + "moving average_wide": self.average_wide, + "initial time_index": self.__init, + } + + def set_params( + self, + carbon_quota: float = None, + demand: Demand = None, + lost: float = None, + optimizer: Optimizer = None, + step: int = 1, + time_index: int = 24 * 7, + time_interval: float = 1, + penalisation: float = None, + carbon_cost: float = None, + plot: str = "default", + average_wide: int = 0, + init: int = 0, + ): + + # init params + self.__params_state = -1 + self.time_index = time_index + self.__init = init + # step is the step of opt budget + self.step = step + self.time_interval = time_interval + self.plot = plot + self.average_wide = average_wide + if demand is not None: + self.set_demand(demand) + if lost is not None: + self.set_constant_lost(lost) + if penalisation is not None: + self.set_penalisation_cost(penalisation) + if carbon_cost is not None: + self.set_carbon_cost(carbon_cost) + if carbon_quota is not None: + self.set_carbon_quota(carbon_quota) + if optimizer is not None: + self.set_optimizer(optimizer) + + def __update_results(self, original, new, init): + output = [] + for step_ in range(0, len(original)): + previous_coef = original[step_]["coef"][:init] + next_coef = new[step_]["coef"] + + previous_loss = self.loss_function(previous_coef, self.time_interval, no_arrange=True) + next_loss = self.loss_function(next_coef, self.time_interval, no_arrange=True) + + production = 0 + u_demand = 0 + carbon_prod = 0 + for i, weighted_coef in enumerate([previous_coef, next_coef]): + for t in range(0, len(weighted_coef)): + production += self.get_production_cost_at_t(weighted_coef[t], t, self.time_interval) + if i == 1: + u_demand += np.abs( + self.get_unsatisfied_demand_at_t(weighted_coef[t], t, self.time_interval, init=init) + ) + else: + u_demand += np.abs(self.get_unsatisfied_demand_at_t(weighted_coef[t], t, self.time_interval)) + carbon_prod += self.get_carbon_production_at_t(weighted_coef[t], self.time_interval) + output.append( + { + "loss": previous_loss + next_loss, + "coef": previous_coef + next_coef, + "production": production, + "unsatisfied demand": u_demand, + "carbon production": carbon_prod, + } + ) + + return output + + ########### + ## PLOTS ## + ########### + def moving_average(self, x, w): + return np.convolve(x, np.ones(w), "valid") / w + + def plotResults(self, optimum: List = [], mode: str = "default", time_interval: float = 1, average_wide: int = 0): + # set the moving average wide + if average_wide == 0: + average_wide = ceil(len(optimum[-1]["coef"]) / 4) + + if mode == "default" or mode == "save" or mode == "best": + # set Y of axe 1 + Y: Dict[str, List[float]] = {} + label_y: List[str] = [] + sum_prod = [] + for c in self.__observable: + label_y.append(c.get_name()) + Y.update({c.get_name(): []}) + if mode != "best": + for array in optimum[-1]["coef"]: + sum_t = 0.0 + for enum, value in enumerate(array): + Y[label_y[enum]].append(value) + sum_t = sum_t + (value * self.get_observable()[enum].get_raw_power()) + sum_prod.append(sum_t) + else: + ##TODO get best loss + pass + + # set Y of axe 2 + Y_: Dict[str, List[float]] = {} + Y_["demand"] = [] + Y_["production"] = sum_prod + for t_i in range(0, len(optimum[-1]["coef"])): + Y_["demand"].append(self.__demand.get_demand_monthly(t_i, time_interval)) + + fig, axs = plt.subplots(1, 2, figsize=(12, 6)) + + # data integration + X = [i for i in range(0, len(optimum[-1]["coef"]))] + for n_axs in range(0, 2): + if n_axs == 0: + for central, value in Y.items(): + smooth_value = self.moving_average(value, average_wide) + axs[n_axs].plot(X[(average_wide - 1) :], smooth_value, ".-", alpha=0.5, lw=2, label=central) + elif n_axs == 1: + axs[n_axs].plot(X, Y_["demand"], label="demand") + smooth_value = self.moving_average(Y_["production"], average_wide) + axs[n_axs].plot(X[(average_wide - 1) :], smooth_value, ".-", alpha=0.5, lw=2, label="production") + + # Add execution_time and loss information + try: + info = ( + "production_cost: " + + "{:.3f}".format(optimum[-1]["loss"]) + + " ; execution_time: " + + "{:.3f}".format(optimum[-1]["elapsed_time"]) + ) + except: + info = "production_cost: " + "{:.3f}".format(optimum[-1]["loss"]) + info += "; demand gap: " + "{:.3f}".format(optimum[-1]["unsatisfied demand"]) + plt.annotate( + info, + xy=(0.5, 0), + xytext=(0, 10), + xycoords=("axes fraction", "figure fraction"), + textcoords="offset points", + size=10, + ha="center", + va="bottom", + ) + + # plots parametrizations + for n_axs in range(0, 2): + axs[n_axs].grid() + axs[n_axs].yaxis.set_tick_params(length=0) + axs[n_axs].xaxis.set_tick_params(length=0) + axs[n_axs].set_xlabel("hours") + # axs[n_axs].yaxis.set_major_formatter(StrMethodFormatter("{x}"+units[0])) + if n_axs == 0: + axs[n_axs].set_ylabel("coef (%)") + elif n_axs == 1: + axs[n_axs].set_ylabel("MWh") + axs[n_axs].legend() + + fig.tight_layout() + try: + path = "results_" + datetime.now().strftime("%d_%m_%Y") + os.makedirs(path) + name = ( + path + + "/" + + "opt_" + + str(self.get_optimizer().get_optimizers()) + + "_" + + datetime.now().strftime("%H%M%S") + + ".png" + ) + fig.savefig(name) + if mode == "default": + plt.show() + + except OSError: + warnings.warn("Can't create the directory " + path + " or already exists") + try: + name = ( + path + + "/" + + "opt_" + + str(self.get_optimizer().get_optimizers()) + + "_" + + datetime.now().strftime("%H%M%S") + + ".png" + ) + fig.savefig(name) + if mode == "default": + plt.show() + except FileNotFoundError: + name = ( + "opt_" + + str(self.get_optimizer().get_optimizers()) + + "_" + + datetime.now().strftime("%H%M%S") + + ".png" + ) + fig.savefig(name) + if mode == "default": + plt.show() + + elif mode == "None": + pass + else: + warnings.warn("Choose an available option : default, save or None") + # plt.show() diff --git a/__init__.py b/mixsimulator/agents/__init__.py similarity index 100% rename from __init__.py rename to mixsimulator/agents/__init__.py diff --git a/mixsimulator/base.py b/mixsimulator/base.py new file mode 100644 index 0000000..486ccd7 --- /dev/null +++ b/mixsimulator/base.py @@ -0,0 +1,23 @@ +import warnings + +from . import Mas_platform as mp +from . import MixSimulator as ms + + +class ElectricityMix: + """ + Regroup all methods to simulate the electricity mix : + - classic + - Multi Agent System (MAS) + """ + + def __init__(self): + pass + + def mix(self, method: str = "classic", carbon_cost: float = 0, penalisation_cost: float = 1000000000000): + if method == "classic": + return ms.MixSimulator(carbon_cost=carbon_cost, penalisation_cost=penalisation_cost) + elif method == "MAS": + return mp.Mas_platform() + else: + warnings.warn("Choose an available option : classic or MAS") diff --git a/mixsimulator/data/.DS_Store b/mixsimulator/data/.DS_Store new file mode 100644 index 0000000..5008ddf Binary files /dev/null and b/mixsimulator/data/.DS_Store differ diff --git a/data/RIToamasina/DIR-TOAMASINA_concat.csv b/mixsimulator/data/RIToamasina/DIR-TOAMASINA_concat.csv similarity index 100% rename from data/RIToamasina/DIR-TOAMASINA_concat.csv rename to mixsimulator/data/RIToamasina/DIR-TOAMASINA_concat.csv diff --git a/data/RIToamasina/dataset_RI_Toamasina.csv b/mixsimulator/data/RIToamasina/dataset_RI_Toamasina.csv similarity index 100% rename from data/RIToamasina/dataset_RI_Toamasina.csv rename to mixsimulator/data/RIToamasina/dataset_RI_Toamasina.csv diff --git a/data/RIToamasina/dataset_RI_Toamasina_v2.csv b/mixsimulator/data/RIToamasina/dataset_RI_Toamasina_v2.csv similarity index 95% rename from data/RIToamasina/dataset_RI_Toamasina_v2.csv rename to mixsimulator/data/RIToamasina/dataset_RI_Toamasina_v2.csv index 77d93f5..95ebfbc 100644 --- a/data/RIToamasina/dataset_RI_Toamasina_v2.csv +++ b/mixsimulator/data/RIToamasina/dataset_RI_Toamasina_v2.csv @@ -1,5 +1,5 @@ centrals;location;tuneable;fuel;hydro;green;fuel_consumption;fuel_cost;availability;init_value;lifetime;carbon_production;raw_power;nb_employees;mean_salary;max_var;height;flow;stock_available;capacity;Demand;lost -Volobe;Volobe;FALSE;FALSE;TRUE;TRUE;1;1;0.650887574;20675460;100;0;6.76;1;1E-10;1;10;10;250;250;11.654;1 +Volobe;Volobe;FALSE;FALSE;TRUE;TRUE;0;0;0.650887574;20675460;100;0;6.76;1;1E-10;1;10;10;250;250;11.654;1 572;Toamasina IV;TRUE;TRUE;FALSE;FALSE;230000;0.000686;0.8333333333;18400;25;730000;6;1;0.000004;0.9;;;;;; 573;Toamasina IV;TRUE;TRUE;FALSE;FALSE;230000;0.000686;0.833333333;18400;25;730000;6;1;0.000004;0.9;;;;;; Tm/ENELEC 1;Betainomby;TRUE;TRUE;FALSE;FALSE;230000;0.000686;0.4938271605;9001.553062985;25;730000;8.1;1;0.000004;0.9;;;;;; diff --git a/mixsimulator/data/RIToamasina/dataset_RI_Toamasina_variation_template.csv b/mixsimulator/data/RIToamasina/dataset_RI_Toamasina_variation_template.csv new file mode 100644 index 0000000..c5cdea7 --- /dev/null +++ b/mixsimulator/data/RIToamasina/dataset_RI_Toamasina_variation_template.csv @@ -0,0 +1,4 @@ +centrals;lower;upper;discrete +572;0.5;1;None +573;0;0.5;0.75:1 +Tm/ENELEC 1;0;0;0:0.25:0.5:0.75:1 diff --git a/data/RIToamasina/demand/toamasina_demand.csv b/mixsimulator/data/RIToamasina/demand/toamasina_demand.csv similarity index 100% rename from data/RIToamasina/demand/toamasina_demand.csv rename to mixsimulator/data/RIToamasina/demand/toamasina_demand.csv diff --git a/data/dataset_CANVAS_fr.xlsx b/mixsimulator/data/dataset_CANVAS_fr.xlsx similarity index 100% rename from data/dataset_CANVAS_fr.xlsx rename to mixsimulator/data/dataset_CANVAS_fr.xlsx diff --git a/data/dataset_demand_CANVAS_fr.xlsx b/mixsimulator/data/dataset_demand_CANVAS_fr.xlsx similarity index 100% rename from data/dataset_demand_CANVAS_fr.xlsx rename to mixsimulator/data/dataset_demand_CANVAS_fr.xlsx diff --git a/data/dataset_variation_CANVAS_fr.xlsx b/mixsimulator/data/dataset_variation_CANVAS_fr.xlsx similarity index 100% rename from data/dataset_variation_CANVAS_fr.xlsx rename to mixsimulator/data/dataset_variation_CANVAS_fr.xlsx diff --git a/centrals/__init__.py b/mixsimulator/demand/__init__.py similarity index 100% rename from centrals/__init__.py rename to mixsimulator/demand/__init__.py diff --git a/mixsimulator/demand/classic/Demand.py b/mixsimulator/demand/classic/Demand.py new file mode 100644 index 0000000..37a4461 --- /dev/null +++ b/mixsimulator/demand/classic/Demand.py @@ -0,0 +1,117 @@ +import csv +import pkgutil +from math import cos, floor, pi +from typing import Any + +import pandas as pd # type: ignore +from prophet import Prophet # type: ignore + + +class Demand: + """ + Manage the Demands data + """ + + def __init__(self, demand: float = 20, var_per_day: float = 0.1, var_per_season: float = 0.1) -> None: + self.__var_per_day = var_per_day + self.__var_per_season = var_per_season + self.__mean_demand = demand + self.__pt = Prophet(seasonality_mode="multiplicative") + self.__periods = 12 * 20 + self.data_demand = None + + def set_forcast_periods(self, periods) -> None: + self.__periods = periods + + def set_data_csv( + self, + bind=None, + raw_data=None, + init_date: str = "2017-01-01", + delimiter: str = ";", + column: str = "Total Ventes", + ): + """ + The method must get a dataset with at least 3 columns + - month : int, + - year : int, + - the monthly demand in kwh (determinated by the parameter "column") + + The method also use a forcast model from prophet to predict future demand. + The periods can be set by set_forcast_periods. + """ + if raw_data is not None: + data = pd.DataFrame(raw_data) + # set columns & index + header = data.iloc[0] + data = data[1:] + data.columns = header + data = data.reset_index(drop=True) + for col in data.columns.tolist(): + try: + # convert numeric values + data[col] = pd.to_numeric(data[col]) + except: + pass + else: + data = pd.read_csv(bind, delimiter=delimiter) + + data["date"] = data["month"].astype("str") + "-" + data["year"].astype("str") + data["datetime"] = pd.to_datetime(data["date"]) + data_to_use = pd.DataFrame() + data_to_use[["datetime", "demands"]] = data.loc[data["datetime"] >= init_date][["datetime", column]] + + # let's forecast + print("############## FORECAST DEMAND ##############") + train = pd.DataFrame() + train[["ds", "y"]] = data[["datetime", column]] + self.__pt.fit(train) + future = self.__pt.make_future_dataframe(periods=self.__periods, freq="MS") + fcst = self.__pt.predict(future) + fcst = fcst.loc[fcst["ds"] > data["datetime"].tail(1).item()] + fcst[["datetime", "demands"]] = fcst[["ds", "yhat"]] + self.data_demand = data_to_use.append(fcst, ignore_index=True) + print("#############################################") + + return self.data_demand + + def set_data_to(self, dataset, delimiter: str = ";"): + data: Any = ... + if dataset == "Toamasina": + # by defaut we keep it "Toamasina" + data = pkgutil.get_data("mixsimulator", "/data/RIToamasina/DIR-TOAMASINA_concat.csv") + data_decoded = csv.reader(data.decode("utf-8").splitlines(), delimiter=delimiter) + self.set_data_csv(raw_data=data_decoded) + else: + # by defaut we keep it "Toamasina" + data = pkgutil.get_data("mixsimulator", "/data/RIToamasina/DIR-TOAMASINA_concat.csv") + data_decoded = csv.reader(data.decode("utf-8").splitlines(), delimiter=delimiter) + self.set_data_csv(raw_data=data_decoded) + + def get_demand(self, t): + self.data_demand.reset_index() + try: + return self.data_demand.iloc[t]["demands"] / 1000 # cause data units in kwh + except IndexError: + print( + "Please check init_date or change the forcasting params : increase periods with Demand.set_forcast_periods" + ) + raise + + def set_mean_demand(self, demand: float): + self.__mean_demand = demand + + def get_demand_approxima(self, t, interval): + demande = self.__mean_demand * ( + 1 + + cos(4 * pi * (t * interval) / 24) * self.__var_per_day + + cos(2 * pi * (t * interval) / (24 * 365)) * self.__var_per_season + ) + return demande + + def get_demand_monthly(self, t, interval): + m = t / (24 * 30) + m = floor(m) + # for now we divide it by 30*24 (better approximation TO DO) + demande = (self.get_demand(m) / (30 * 24)) * (1 + cos(4 * pi * (t * interval) / 24) * self.__var_per_day) + return demande diff --git a/nevergradBased/__init__.py b/mixsimulator/demand/classic/__init__.py similarity index 100% rename from nevergradBased/__init__.py rename to mixsimulator/demand/classic/__init__.py diff --git a/mixsimulator/demand/mas/Demand.py b/mixsimulator/demand/mas/Demand.py new file mode 100644 index 0000000..c75e243 --- /dev/null +++ b/mixsimulator/demand/mas/Demand.py @@ -0,0 +1,150 @@ +import csv +import json +import pkgutil +from math import cos, floor, pi +from typing import Any, List + +import pandas as pd # type: ignore +from prophet import Prophet # type: ignore + +from ...agents.Agent import Agent + + +def massive_divergence(demande, demande2): + return True + + +class Demand(Agent): + def __init__(self, model=None, demand: float = 20, var_per_day: float = 0.1, var_per_season: float = 0.1) -> None: + super().__init__() + if model is not None: + self.set_model(model) + self.set_type("Demand") + self.__demand: List[float] = [] + + ### from old demand implementation + self.__var_per_day = var_per_day + self.__var_per_season = var_per_season + self.__mean_demand = demand + self.__pt = Prophet(seasonality_mode="multiplicative") + self.__periods = 12 * 20 + self.data_demand = None + + ### COMMUNICATION + def __notify_demand_value_change(self): + signal = json.loads(open(self._code_files))["8080"] + signal["values"] = self.__demand + signal["t_from"] = 3 + self._notify_observers(signal) + + def set_demande_change(self, demand: List[float]) -> None: + errored_demand = self.__demand + self.__demand = demand + if massive_divergence(errored_demand, self.__demand): + self.__notify_demand_value_change() + + ### PARAMETRIZATION + def set_model(self, model) -> None: + self.__model = model + + def predict_demand(self, time_series: pd.DataFrame) -> None: + self.__demand = list(self.__model.predict(time_series)) + + ### TMP WITH PROPHET + def set_forecast_periods(self, periods) -> None: + self.__periods = periods + + def forecast_with_prophet( + self, + bind=None, + raw_data=None, + init_date: str = "2017-01-01", + delimiter: str = ";", + column: str = "Total Ventes", + ): + """ + The method must get a dataset with at least 3 columns + - month : int, + - year : int, + - the monthly demand in kwh (determinated by the parameter "column") + + The method also use a forcast model from prophet to predict future demand. + The periods can be set by set_forcast_periods. + """ + if raw_data is not None: + data = pd.DataFrame(raw_data) + # set columns & index + header = data.iloc[0] + data = data[1:] + data.columns = header + data = data.reset_index(drop=True) + for col in data.columns.tolist(): + try: + # convert numeric values + data[col] = pd.to_numeric(data[col]) + except: + pass + else: + data = pd.read_csv(bind, delimiter=delimiter) + + data["date"] = data["month"].astype("str") + "-" + data["year"].astype("str") + data["datetime"] = pd.to_datetime(data["date"]) + data_to_use = pd.DataFrame() + data_to_use[["datetime", "demands"]] = data.loc[data["datetime"] >= init_date][["datetime", column]] + + # let's forecast + print("############## FORECAST DEMAND ##############") + train = pd.DataFrame() + train[["ds", "y"]] = data[["datetime", column]] + self.__pt.fit(train) + future = self.__pt.make_future_dataframe(periods=self.__periods, freq="MS") + fcst = self.__pt.predict(future) + fcst = fcst.loc[fcst["ds"] > data["datetime"].tail(1).item()] + fcst[["datetime", "demands"]] = fcst[["ds", "yhat"]] + self.data_demand = data_to_use.append(fcst, ignore_index=True) + print("#############################################") + + return self.data_demand + + def set_data_to(self, dataset, delimiter: str = ";"): + data: Any = ... + if dataset == "Toamasina": + # by defaut we keep it "Toamasina" + data = pkgutil.get_data("mixsimulator", "/data/RIToamasina/DIR-TOAMASINA_concat.csv") + data_decoded = csv.reader(data.decode("utf-8").splitlines(), delimiter=delimiter) + self.forecast_with_prophet(raw_data=data_decoded) + else: + # by defaut we keep it "Toamasina" + data = pkgutil.get_data("mixsimulator", "/data/RIToamasina/DIR-TOAMASINA_concat.csv") + data_decoded = csv.reader(data.decode("utf-8").splitlines(), delimiter=delimiter) + self.forecast_with_prophet(raw_data=data_decoded) + + def get_demand(self, t): + self.data_demand.reset_index() + try: + return self.data_demand.iloc[t]["demands"] / 1000 # cause data units in kwh + except IndexError: + print( + "Please check init_date or change the forcasting params : increase periods with Demand.set_forcast_periods" + ) + raise + + def set_mean_demand(self, demand: float): + self.__mean_demand = demand + + def get_demand_approxima(self, t, interval): + demande = self.__mean_demand * ( + 1 + + cos(4 * pi * (t * interval) / 24) * self.__var_per_day + + cos(2 * pi * (t * interval) / (24 * 365)) * self.__var_per_season + ) + return demande + + def get_demand_monthly(self, t, interval, init: int = 0): + m = (t + init) / (24 * 30) + m = floor(m) + # for now we divide it by 30*24 (better approximation TO DO) + demande = (self.get_demand(m) / (30 * 24)) * ( + 1 + cos(4 * pi * ((t + init) * interval) / 24) * self.__var_per_day + ) + return demande diff --git a/mixsimulator/demand/mas/__init__.py b/mixsimulator/demand/mas/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/mixsimulator/nevergradBased/Optimizer.py b/mixsimulator/nevergradBased/Optimizer.py new file mode 100644 index 0000000..8276b76 --- /dev/null +++ b/mixsimulator/nevergradBased/Optimizer.py @@ -0,0 +1,158 @@ +import time +import tqdm +from typing import List, Type + +import nevergrad as ng +import numpy as np # type: ignore + + +class Optimizer: + """ + adaptation of nevergrad optimizers to the project + auto-parametrization + + list of available: self.getOptimizerList() + + """ + + def __init__( + self, + opt=[ng.optimizers.OnePlusOne], + budget: List[int] = [100], + num_worker: int = 1, + instrum=ng.p.Array(shape=(2,)), + ): + self.set_budget(budget) + self.set_parametrization(instrum) + self.set_num_worker(num_worker) + + ### available optimizers + self.__available_optimizers = {Type[str]: Type[object]} + nevergrad_optimizers = list(ng.optimizers.registry.keys()) + for ng_opt in nevergrad_optimizers: + self.__available_optimizers.update({ng_opt: ng.optimizers.registry[ng_opt]}) + + convert_opt = [] + for opt_ng in opt: + if opt_ng in list(self.__available_optimizers.values()): + convert_opt.append(opt_ng) + continue + try: + convert_opt.append(self.__available_optimizers[str(opt_ng)]) + except: + raise KeyError(str(opt_ng) + " not available") + self.__optimizers = convert_opt + + def get_available_optimizers(self): + tmp = Optimizer() + return tmp.__available_optimizers.keys() + + def get_unavailable_optimizers(self): + tmp = Optimizer() + result = [] + tmp = tmp.__available_optimizers.keys() + list = sorted(ng.optimizers.registry.keys()) + for opt in list: + if opt not in tmp: + result.append(opt) + return result + + def get_optimizers(self): + return self.__optimizers + + def set_parametrization(self, instrum): + self.__parametrization = instrum + + def get_parametrization(self): + return self.__parametrization + + def set_budget(self, budget: List[int] = [100]): + # setting budget + # possible to do? auto calculate a optimal budget depending of the size of data + self.__budget = budget + + def get_budget(self): + return self.__budget + + def set_num_worker(self, n_work): + # setting budget + # possible to do? auto calculate a optimal budget depending of the size of data + self.__num_worker = n_work + + def get_num_worker(self): + return self.__num_worker + + def get_params(self) -> dict: + return { + "optimizer(s)": self.__optimizers, + "budget(s)": self.get_budget(), + "parametrization": self.get_parametrization(), + "num_worker": self.get_num_worker(), + } + + def optimize(self, mix=None, func_to_optimize=None, constraints=None, step: int = 1, init: int = 0): + + # setting budgets + budgets = self.get_budget() + + if len(budgets) != len(self.__optimizers): + raise IndexError("\n The length of budgets and the length of optimizers list should be the same.\n") + + total_budget = 0 + for algo_budget in budgets: + total_budget += algo_budget + result = [] + start_time = time.time() + + # optimization under constraints + chaining_algo = ng.optimizers.Chaining(self.__optimizers, budgets[:-1]) + optimizer = chaining_algo( + parametrization=self.get_parametrization(), budget=budgets[-1], num_workers=self.get_num_worker() + ) + + # let's minimize + print("Optimization over {} budgets :".format(total_budget)) + for tmp_budget in tqdm.tqdm(range(0, total_budget)): + x = optimizer.ask() + loss = func_to_optimize(*x.args, constraints["time_interval"], init=init) + optimizer.tell(x, loss) + if (tmp_budget + 1) % step == 0: + result_per_budget = {} + recommendation = optimizer.provide_recommendation() + result_per_budget.update( + {"loss": func_to_optimize(recommendation.value, constraints["time_interval"])} + ) + result_per_budget.update({"coef": recommendation.value}) + + # Readjustment + if mix is not None: + usage_coef = mix.arrange_coef_as_array_of_array(recommendation.value) + try: + weighted_coef = mix.get_weighted_coef( + usage_coef, time_interval=constraints["time_interval"], init=init + ) + except: + weighted_coef = mix.get_weighted_coef(usage_coef, time_interval=constraints["time_interval"]) + production = 0 + u_demand = 0 + carbon_prod = 0 + for t in range(0, len(weighted_coef)): + production += mix.get_production_cost_at_t(weighted_coef[t], t, constraints["time_interval"]) + try: + u_demand += np.abs( + mix.get_unsatisfied_demand_at_t( + weighted_coef[t], t, constraints["time_interval"], init=init + ) + ) + except: + u_demand += np.abs( + mix.get_unsatisfied_demand_at_t(weighted_coef[t], t, constraints["time_interval"]) + ) + carbon_prod += mix.get_carbon_production_at_t(weighted_coef[t], constraints["time_interval"]) + result_per_budget.update({"production": production}) + result_per_budget.update({"unsatisfied demand": u_demand}) + result_per_budget.update({"carbon production": carbon_prod}) + result_per_budget.update({"elapsed_time": time.time() - start_time}) + result.append(result_per_budget) + + # TODO --> Check Constraints + return result diff --git a/mixsimulator/nevergradBased/__init__.py b/mixsimulator/nevergradBased/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/mixsimulator/params_files/exchange_code.json b/mixsimulator/params_files/exchange_code.json new file mode 100644 index 0000000..df1e7af --- /dev/null +++ b/mixsimulator/params_files/exchange_code.json @@ -0,0 +1,36 @@ +{ + "100" : { + "type" : "register", + "id" : null, + "t_from": null, + "code": 100 + }, + "200" : { + "status" : "OK", + "type" : "health", + "id" : null, + "t_from": null, + "code" : 200 + }, + "400" : { + "status" : "DOWN", + "type" : "health", + "id" : null, + "t_from": null, + "code" : 400 + }, + + "101" : { + "type" : "disponibility", + "value" : null, + "id" : null, + "t_from": null, + "code" : 101 + }, + + "8080" : { + "type" : "demand_change", + "t_from" : null, + "values" : null + } +} \ No newline at end of file diff --git a/mixsimulator/params_files/settings.json b/mixsimulator/params_files/settings.json new file mode 100644 index 0000000..fdb14dc --- /dev/null +++ b/mixsimulator/params_files/settings.json @@ -0,0 +1,12 @@ +{ + "power_plant_environement": + { + "simple_weather": { + "url": "https://api.openweathermap.org/data/2.5/weather?lat={lat}&lon={lon}&appid=3159336af6f662a719a184f1bb722347" + }, + "precipitation": { + "url": "https://maps.openweathermap.org/maps/2.0/weather/1h/PR0/{zoom}/{x_tile}/{y_tile}?appid=3159336af6f662a719a184f1bb722347" + } + } + +} \ No newline at end of file diff --git a/mixsimulator/power_plants/__init__.py b/mixsimulator/power_plants/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/centrals/HydroCentral.py b/mixsimulator/power_plants/classic/HydroCentral.py similarity index 60% rename from centrals/HydroCentral.py rename to mixsimulator/power_plants/classic/HydroCentral.py index b15ad29..6c9960a 100644 --- a/centrals/HydroCentral.py +++ b/mixsimulator/power_plants/classic/HydroCentral.py @@ -1,18 +1,28 @@ -from .PowerCentral import PowerCentral as pc from math import cos, pi -from typing import List + +from .PowerCentral import PowerCentral as pc + class HydroCentral(pc): """ - Class of power plant with the specifications of a HydroElectric Power Plant + Class of power plant with the specifications of a HydroElectric Power Plant """ - def __init__(self, hauteur: float, moyenne_apport: float, capacity: float, available_stock: float, var_per_day: float, var_per_season: float ) -> None: + + def __init__( + self, + hauteur: float, + moyenne_apport: float, + capacity: float, + available_stock: float, + var_per_day: float, + var_per_season: float, + ) -> None: super().__init__() self.__hauteur = hauteur - self.__moyenne_apport = moyenne_apport # m3/s - self.__capacity = capacity # m3 + self.__moyenne_apport = moyenne_apport # m3/s + self.__capacity = capacity # m3 self.__init_stock = available_stock - self.__available_stock = [available_stock]# m3 + self.__available_stock = [available_stock] # m3 self.__tuneable = True self.__var_per_day = var_per_day self.__var_per_season = var_per_season @@ -21,14 +31,18 @@ def isHydro(self) -> bool: return True def __get_natural_availability(self, t) -> float: - debit_t = self.__moyenne_apport * (1 + (cos(2 * pi * ( t )/ 24))* self.__var_per_day + (cos(2 * pi * ( t ) / (24*365)))* self.__var_per_season) - power = (self.__hauteur * debit_t * 9.8 * 0.9 * 997)/1000000 # unit is MW - return power/self._raw_power + debit_t = self.__moyenne_apport * ( + 1 + + (cos(2 * pi * (t) / 24)) * self.__var_per_day + + (cos(2 * pi * (t) / (24 * 365))) * self.__var_per_season + ) + power = (self.__hauteur * debit_t * 9.8 * 0.9 * 997) / 1000000 # unit is MW + return power / self._raw_power - def __get__artificial_availability(self, t) ->float: - debit_t_max = self.__get_available_stock(t)/(3600) - power = (self.__hauteur * debit_t_max * 9.8 * 0.9 * 997) /1000000 # unit is MW - return power/self._raw_power + def __get__artificial_availability(self, t) -> float: + debit_t_max = self.__get_available_stock(t) / (3600) + power = (self.__hauteur * debit_t_max * 9.8 * 0.9 * 997) / 1000000 # unit is MW + return power / self._raw_power def get_availability(self, t) -> float: dummy_availability = self.__get__artificial_availability(t) + self.__get_natural_availability(t) @@ -44,7 +58,7 @@ def back_propagate(self, usage_coef, t, interval): diff = usage_coef - self.__get_natural_availability(t) diff_power = diff * self._raw_power # back to m3/s so diff_power has to be in W not in MW ===> *1000 - bandwidth = (diff_power* 1000000) / (9.8*self.__hauteur* 0.9 * 997) + bandwidth = (diff_power * 1000000) / (9.8 * self.__hauteur * 0.9 * 997) self.__update_stock(bandwidth, interval, t) def __update_stock(self, bandwidth, interval, t): @@ -52,7 +66,7 @@ def __update_stock(self, bandwidth, interval, t): #### if positiv, we have used stocked water #### if not the stock has been increased current_availability = self.__available_stock[t] - current_availability += -bandwidth*(interval*3600) + current_availability += -bandwidth * (interval * 3600) if current_availability > self.__capacity: current_availability = self.__capacity self.__available_stock.append(current_availability) @@ -62,4 +76,4 @@ def reset_central(self): self.__available_stock = [self.__init_stock] def get_stock_evolution(self): - return self.__available_stock \ No newline at end of file + return self.__available_stock diff --git a/centrals/PowerCentral.py b/mixsimulator/power_plants/classic/PowerCentral.py similarity index 60% rename from centrals/PowerCentral.py rename to mixsimulator/power_plants/classic/PowerCentral.py index ff8ebec..2e3539c 100644 --- a/centrals/PowerCentral.py +++ b/mixsimulator/power_plants/classic/PowerCentral.py @@ -1,29 +1,32 @@ -from typing import List +from typing import Any, List + import nevergrad as ng + class PowerCentral: """ - Class for basic power plant - it has all the common parameters of the control units (central) + Class for basic power plant + it has all the common parameters of the control units (central) """ - def __init__(self, tuneable:bool=False) -> None: - self._id = "0" - self._changeRate = 0. #(percent) - self._initial_value = 0. - self._lifetime = 0 #in hour - self._carbon_prod = 0. #g/MWh - self._raw_power = 0. #MW - self._availability = 1. #of the source + + def __init__(self, tuneable: bool = False) -> None: + self._id = "0" + self._changeRate = 0.0 # (percent) + self._initial_value = 0.0 + self._lifetime = 0 # in hour + self._carbon_prod = 0.0 # g/MWh + self._raw_power = 0.0 # MW + self._availability = 1.0 # of the source self._nb_employes = 1 - self._mean_salary = 0. #per month + self._mean_salary = 0.0 # per month self._tuneable = tuneable - self.__fuel_cost = 0. #$/g - self.__fuel_consumption = 0. #g/MWh + self.__fuel_cost = 0.0 # $/g + self.__fuel_consumption = 0.0 # g/MWh self._init_cur_usage = 0 self._cur_usage = 0 self._max_var = 1 - self._lower = 0. - self._upper = 1. + self._lower: Any = 0.0 + self._upper: Any = 1.0 self._choices = None def set_init_cur_usage(self, init_usage): @@ -64,25 +67,25 @@ def get_amortized_cost(self, time_index): if time_index > self._lifetime * 365 * 24: return 0 else: - return ( self._initial_value / ( self._lifetime * 365 * 24 ) ) + return self._initial_value / (self._lifetime * 365 * 24) def is_tuneable(self) -> bool: - #is controlled or not + # is controlled or not return self._tuneable - def get_carbon_production(self) -> float: # g/MWh + def get_carbon_production(self) -> float: # g/MWh return self._carbon_prod - def set_carbon_prod(self, carbonCost: float=0) -> None: + def set_carbon_prod(self, carbonCost: float = 0) -> None: self._carbon_prod = carbonCost def set_raw_power(self, rawPower): self._raw_power = rawPower - def get_raw_power(self) -> float: # MW + def get_raw_power(self) -> float: # MW return self._raw_power - def get_availability(self, time_index) -> float: # percent + def get_availability(self, time_index) -> float: # percent return self._availability def reset_central(self): @@ -91,24 +94,24 @@ def reset_central(self): def get_max_availability(self, time_index) -> float: theorical_availability = self.get_availability(time_index) if self._cur_usage + self._max_var <= theorical_availability: - theorical_availability = self._cur_usage + self._max_var + theorical_availability = self._cur_usage + self._max_var return theorical_availability def get_min_availability(self, time_index) -> float: theorical_availability = 0 if self._cur_usage - self._max_var >= theorical_availability: - theorical_availability = self._cur_usage - self._max_var + theorical_availability = self._cur_usage - self._max_var return theorical_availability def back_propagate(self, usage_coef, t, time_interval): self._cur_usage = usage_coef - + def set_availability(self, availability: float): self._availability = availability - def get_employees_salary(self, total_working_time_per_day=8) ->float: - #mean salary per hour - perHourMeanSalary = self._mean_salary/(31*total_working_time_per_day) + def get_employees_salary(self, total_working_time_per_day=8) -> float: + # mean salary per hour + perHourMeanSalary = self._mean_salary / (31 * total_working_time_per_day) return perHourMeanSalary * self._nb_employes def set_mean_employees_salary(self, mean_salary): @@ -118,38 +121,42 @@ def set_tuneable(self, tuneable: bool) -> None: self._tuneable = tuneable def __getUsageCoef(self, usage_coef: float) -> None: - if(self._tuneable): + if self._tuneable: usage_coef = min(self._availability, usage_coef) else: usage_coef = self._availability - - def set_variation_params(self, lower: float, upper : float, choices : List[float] = None) -> None: + + def set_variation_params(self, lower: Any, upper: Any, choices: Any = None) -> None: self._lower = lower self._upper = upper self._choices = choices - + def get_variation_params(self) -> ng.p.Choice: - final_params = [] - if self._choices != None : - if self._lower == self._upper : + final_params = [] + if self._choices is not None: + if self._lower == self._upper: return ng.p.Choice(self._choices) - else : - for low, up in zip(self._lower, self._upper): - if low == up : continue - final_params.append(ng.p.Scalar(lower=low ,upper=up)) + else: + lows: Any = [self._lower] + uppers: Any = [self._upper] + for low, up in zip(lows, uppers): + if low[0] == up[0]: + continue + final_params.append(ng.p.Scalar(lower=low[0], upper=up[0])) discret = ng.p.Choice(self._choices) final_params.append(discret) params = ng.p.Choice(final_params) return params - - else : - if self._lower != 0. and self._upper != 1. : - for low, up in zip(self._lower, self._upper): - final_params.append(ng.p.Scalar(lower=low ,upper=up)) + + else: + if self._lower != 0.0 and self._upper != 1.0: + lows = [self._lower] + uppers = [self._upper] + for low, up in zip(lows, uppers): + final_params.append(ng.p.Scalar(lower=low[0], upper=up[0])) params = ng.p.Choice(final_params) - return params - else : - final_params.append(ng.p.Scalar(lower=self._lower ,upper=self._upper)) + return params + else: + final_params.append(ng.p.Scalar(lower=self._lower, upper=self._upper)) params = ng.p.Choice(final_params) return params - \ No newline at end of file diff --git a/mixsimulator/power_plants/classic/__init__.py b/mixsimulator/power_plants/classic/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/mixsimulator/power_plants/mas/Hydropowerplant.py b/mixsimulator/power_plants/mas/Hydropowerplant.py new file mode 100644 index 0000000..79e80a7 --- /dev/null +++ b/mixsimulator/power_plants/mas/Hydropowerplant.py @@ -0,0 +1,89 @@ +from math import cos, pi + +from .PowerPlant import PowerPlant + + +class Hydropowerplant(PowerPlant): + """ + Agent simulating a Hydroelectic facility + """ + + def __init__( + self, + hauteur: float, + moyenne_apport: float, + capacity: float, + available_stock: float, + var_per_day: float, + var_per_season: float, + ) -> None: + super().__init__() + self.set_type("Hydropowerplant") + self.__hauteur = hauteur + self.__moyenne_apport = moyenne_apport # m3/s + self.__capacity = capacity # m3 + self.__init_stock = available_stock + self.__available_stock = [available_stock] # m3 + self.__tuneable = True + self.__var_per_day = var_per_day + self.__var_per_season = var_per_season + + def isHydro(self) -> bool: + return True + + def __get_natural_availability(self, t, init: int = 0) -> float: + debit_t = self.__moyenne_apport * ( + 1 + + (cos(2 * pi * (t + init) / 24)) * self.__var_per_day + + (cos(2 * pi * (t + init) / (24 * 365))) * self.__var_per_season + ) + power = (self.__hauteur * debit_t * 9.8 * 0.9 * 997) / 1000000 # unit is MW + if self._raw_power == 0.0: + return 0.0 + else: + return power / self._raw_power + + def __get__artificial_availability(self, t, init: int = 0) -> float: + debit_t_max = self.__get_available_stock(t + init) / (3600) + power = (self.__hauteur * debit_t_max * 9.8 * 0.9 * 997) / 1000000 # unit is MW + if self._raw_power == 0.0: + return 0.0 + else: + return power / self._raw_power + + def get_availability(self, t, init: int = 0) -> float: + dummy_availability = self.__get__artificial_availability(t, init=init) + self.__get_natural_availability( + t, init=init + ) + if dummy_availability > 1: + dummy_availability = 1 + return dummy_availability + + def __get_available_stock(self, t): + return self.__available_stock[t] + + def back_propagate(self, usage_coef, t, interval, init: int = 0): + super().back_propagate(usage_coef, t, interval, init=init) + diff = usage_coef - self.__get_natural_availability(t, init=init) + diff_power = diff * self._raw_power + # back to m3/s so diff_power has to be in W not in MW ===> *1000 + bandwidth = (diff_power * 1000000) / (9.8 * self.__hauteur * 0.9 * 997) + self.__update_stock(bandwidth, interval, t, init=init) + + def __update_stock(self, bandwidth, interval, t, init: int = 0): + # bandwidth can be either positive or negative + #### if positiv, we have used stocked water + #### if not the stock has been increased + current_availability = self.__available_stock[t + init] + current_availability += -bandwidth * (interval * 3600) + if current_availability > self.__capacity: + current_availability = self.__capacity + if init == 0: + self.__available_stock.append(current_availability) + else: + self.__available_stock[t + init] = current_availability + + def get_stock_evolution(self): + return self.__available_stock + + ### TODO : add/improve prediction function, data monitoring and signal system diff --git a/mixsimulator/power_plants/mas/PowerPlant.py b/mixsimulator/power_plants/mas/PowerPlant.py new file mode 100644 index 0000000..22437f6 --- /dev/null +++ b/mixsimulator/power_plants/mas/PowerPlant.py @@ -0,0 +1,223 @@ +# import sklearn.linear_model as linear_model +import json +import math +import pkgutil +from typing import Any, Dict, List + +import nevergrad as ng + +from ...agents.Agent import Agent + + +class PowerPlant(Agent): + """ + Agent simulating all the common functions of a power plant + """ + + def __init__(self, tuneable: bool = False) -> None: + super().__init__() + data_json: Any = pkgutil.get_data("mixsimulator", "params_files/settings.json") + self.__api_setting = json.loads(data_json.decode("utf-8")) + self.__changeRate = 0.0 # (percent) + self.__initial_value = 0.0 + self.__lifetime = 0 # in hour + self.__carbon_prod = 0.0 # g/MWh + self._raw_power = 0.0 # MW + self._original_raw_power = 0.0 + self.__availability = 1.0 # of the source + self.__nb_employes = 1 + self.__mean_salary = 0.0 # per month + self.__tuneable = tuneable + self.__fuel_cost = 0.0 # $/g + self.__fuel_consumption = 0.0 # g/MWh + self.__init_cur_usage = 0 + self.__cur_usage = 0 + self.__max_var = 1 + self._lower: Any = 0.0 + self._upper: Any = 1.0 + self._choices = None + self.__lat: float = 0.0 + self.__long: float = 0.0 + self.__x_tile: int = 0 + self.__y_tile: int = 0 + self.__zoom: float = 0.0 + + ### ENVIRONEMENT MONITORING & MODELISATION + def set_location(self, location: Dict, zoom: int = 2) -> None: + self.__lat = location["lat"] + self.__long = location["long"] + lat_rad: Any = math.radians(self.__lat) + n = 2.0**zoom + self.__x_tile = int((self.__long + 180.0) / 360.0 * n) + self.__y_tile = int((1.0 - math.asinh(math.tan(lat_rad)) / math.pi) / 2.0 * n) + self.__zoom = zoom + self._schedule_action({self.get_environement: 60}) + + def get_location(self) -> Dict: + return { + "lat": self.__lat, + "long": self.__long, + "x_tile": self.__x_tile, + "y_tile": self.__y_tile, + "zoom": self.__zoom, + } + + def get_environement(self) -> str: + url = self.__api_setting["power_plant_environement"]["simple_weather"]["url"] + ### TO DO + # FETCH FROM API + print("getting env") + + return "value" + + ### COMMUNICATION + def _notify_is_up(self, t_from=0) -> None: + data_json: Any = pkgutil.get_data("mixsimulator", self._code_files) + signal = json.loads(data_json.decode("utf-8"))["200"] + signal["id"] = self.get_id() + signal["t_from"] = t_from + self._notify_observers(signal) + + def _notify_is_down(self, t_from=0) -> None: + data_json: Any = pkgutil.get_data("mixsimulator", self._code_files) + signal = json.loads(data_json.decode("utf-8"))["400"] + signal["id"] = self.get_id() + signal["t_from"] = t_from + self._notify_observers(signal) + + def _notify_disponibility(self) -> None: + pass + + ### PARAMETRIZATION + def set_init_cur_usage(self, init_usage): + self.__init_cur_usage = init_usage + self.__cur_usage = init_usage + + def set_max_var(self, max_var): + self.__max_var = max_var + + def set_nb_employees(self, nb_employees): + self.__nb_employes = nb_employees + + def set_initial_value(self, initial_value): + self.__initial_value = initial_value + + def set_lifetime(self, lifetime): + self.__lifetime = lifetime + + def set_fuel_consumption(self, fuel_consumption): + self.__fuel_consumption = fuel_consumption + + def get_fuel_consumption(self): + return self.__fuel_consumption + + def set_fuel_cost(self, fuel_cost): + self.__fuel_cost = fuel_cost + + def get_fuel_cost(self): + return self.__fuel_cost + + def get_amortized_cost(self, time_index): + if time_index > self.__lifetime * 365 * 24: + return 0 + else: + return self.__initial_value / (self.__lifetime * 365 * 24) + + def is_tuneable(self) -> bool: + # is controlled or not + return self.__tuneable + + def get_carbon_production(self) -> float: # g/MWh + return self.__carbon_prod + + def set_carbon_prod(self, carbonCost: float = 0) -> None: + self.__carbon_prod = carbonCost + + def set_raw_power(self, rawPower): + self._raw_power = rawPower + + def get_raw_power(self) -> float: # MW + return self._raw_power + + def get_availability(self, time_index, init: int = 0) -> float: # percent + return self.__availability + + def reset_powerplant(self): + self.__cur_usage = self.__init_cur_usage + + def get_max_availability(self, time_index, init: int = 0) -> float: + theorical_availability = self.get_availability(time_index, init=init) + if self.__cur_usage + self.__max_var <= theorical_availability: + theorical_availability = self.__cur_usage + self.__max_var + return theorical_availability + + def get_min_availability(self) -> float: + theorical_availability = 0 + if self.__cur_usage - self.__max_var >= theorical_availability: + theorical_availability = self.__cur_usage - self.__max_var + return theorical_availability + + def back_propagate(self, usage_coef, t, time_interval, init: int = 0): + self.__cur_usage = usage_coef + + def set_availability(self, availability: float): + self.__availability = availability + + def get_employees_salary(self, total_working_time_per_day=8) -> float: + # mean salary per hour + perHourMeanSalary = self.__mean_salary / (31 * total_working_time_per_day) + return perHourMeanSalary * self.__nb_employes + + def set_mean_employees_salary(self, mean_salary): + self.__mean_salary = mean_salary + + def set_tuneable(self, tuneable: bool) -> None: + self.__tuneable = tuneable + + def shutdown(self) -> None: + self._original_raw_power = self._raw_power + self._raw_power = 0.0 + + def power_on(self) -> None: + self._raw_power = self._original_raw_power + + def __getUsageCoef(self, usage_coef: float) -> None: + if self.__tuneable: + usage_coef = min(self.__availability, usage_coef) + else: + usage_coef = self.__availability + + def set_variation_params(self, lower: Any, upper: Any, choices: Any = None) -> None: + self._lower = lower + self._upper = upper + self._choices = choices + + def get_variation_params(self) -> ng.p.Choice: + final_params: List[Any] = [] + if self._choices is not None: + if self._lower == self._upper: + return ng.p.Choice(self._choices) + else: + lows: Any = [self._lower] + uppers: Any = [self._upper] + for low, up in zip(lows, uppers): + if low == up: + continue + final_params.append(ng.p.Scalar(lower=low, upper=up)) + discret = ng.p.Choice(self._choices) + final_params.append(discret) + params = ng.p.Choice(final_params) + return params + + else: + if self._lower != 0.0 and self._upper != 1.0: + lows = [self._lower] + uppers = [self._upper] + for low, up in zip(lows, uppers): + final_params.append(ng.p.Scalar(lower=low, upper=up)) + params = ng.p.Choice(final_params) + return params + else: + final_params.append(ng.p.Scalar(lower=self._lower, upper=self._upper)) + params = ng.p.Choice(final_params) + return params diff --git a/mixsimulator/power_plants/mas/Thermalpowerplant.py b/mixsimulator/power_plants/mas/Thermalpowerplant.py new file mode 100644 index 0000000..7e3e994 --- /dev/null +++ b/mixsimulator/power_plants/mas/Thermalpowerplant.py @@ -0,0 +1,16 @@ +from .PowerPlant import PowerPlant + + +class Thermalpowerplant(PowerPlant): + """ + Agent simulating a fuel, nuclear or solar thermal power plant + """ + + def __init__(self): + super().__init__() + self.set_type("Thermalpowerplant") + + def __get_water_availability(self): + pass + + ### TODO : add/improve prediction function, data monitoring and signal system diff --git a/mixsimulator/power_plants/mas/__init__.py b/mixsimulator/power_plants/mas/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/test.py b/mixsimulator/test_classic.py old mode 100755 new mode 100644 similarity index 61% rename from test.py rename to mixsimulator/test_classic.py index db9889b..1e43d19 --- a/test.py +++ b/mixsimulator/test_classic.py @@ -1,43 +1,46 @@ -import mixsimulator.MixSimulator as ms -import mixsimulator.nevergradBased.Optimizer as opt +import mixsimulator.MixSimulator as ms +import mixsimulator.nevergradBased.Optimizer as opt from mixsimulator.Evaluation import EvaluationBudget -import mixsimulator.Demand as dm -""" -(1) Configure nevergrad optimizers +from .demand.classic.Demand import Demand + +""" +(1) Configure nevergrad optimizers Default Parameters : ---------- - opt = [ng.optimizers.OnePlusOne], - budget: List[int] = [100], - num_worker: int = 1, + opt = [ng.optimizers.OnePlusOne], + budget: List[int] = [100], + num_worker: int = 1, instrum = ng.p.Array(shape=(2,)) """ -opt_CMA = opt.Optimizer(opt = ["CMA"], budget = [20], num_worker = 1) -opt_CMA_30 = opt.Optimizer(opt = ["CMA"], budget = [20], num_worker = 30) +opt_CMA = opt.Optimizer(opt=["CMA"], budget=[20], num_worker=1) +opt_CMA_30 = opt.Optimizer(opt=["CMA"], budget=[20], num_worker=30) -""" -(2) Init MixSimulator instance """ -mix = ms.MixSimulator() +(2) Init MixSimulator instance +""" +mix = ms() """ (3) Set dataset to use - - all datasets available on : + - all datasets available on : https://github.com/Foloso/MixSimulator/tree/master/data/RIToamasina """ -mix.set_data_csv("data/RIToamasina/dataset_RI_Toamasina_v2.csv",delimiter=";") +# mix.set_data_csv("data/RIToamasina/dataset_RI_Toamasina_v2.csv", delimiter=";") -### or use it for default dataset (RI_Toamasina version 2) -#mix.set_data_to("Toamasina") +""" +or use it for default dataset (RI_Toamasina version 2) +""" +mix.set_data_to("Toamasina") """ -(4) For variation limits dataset, there is not yet default dataset +(4) For variation limits dataset, there is not yet default dataset - the "dataset_RI_Toamasina_variation_template.csv" is a random dataset - - Description (fr): + - Description (fr): COLONNES DESCRIPTION -------- -------------------------------------------------------- centrals Nom ou identifiant de la centrale @@ -46,38 +49,38 @@ discrete Les valeurs discrètes des variations fixes des centrales """ -mix.set_variation_csv("data/RIToamasina/dataset_RI_Toamasina_variation_template.csv",delimiter=";") +# mix.set_variation_csv("data/RIToamasina/dataset_RI_Toamasina_variation_template.csv", delimiter=";") """ (5) Load the dataset demand to use : based on a monthly data (column "Total Ventes") from 2008 to 2017 Beyond that date it's a forecasting by prophet! - - data available on : + - data available on : https://github.com/Foloso/MixSimulator/tree/master/data/RIToamasina/demand """ -demand = dm.Demand() -data_demand = demand.set_data_csv("data/RIToamasina/DIR-TOAMASINA_concat.csv", delimiter = ",") +demand = Demand() +# data_demand = demand.set_data_csv("data/RIToamasina/DIR-TOAMASINA_concat.csv", delimiter=",") +### or for default dataset +demand.set_data_to("Toamasina", delimiter=",") """ The method must get a dataset with at least 3 columns - - month : int, + - month : int, - year : int, - the monthly demand in kwh (determinated by the parameter "column") - + The method also use a forcast model from prophet to predict future demand. The periods can be set by set_forcast_periods. """ -#or for default dataset -#demand.set_data_to("Toamasina",delimiter=",") -""" -(6) Then set the demand +""" +(6) Then set the demand """ mix.set_demand(demand) -""" +""" (7) OPTIMIZATION : - + Default parameters : --------------------------- carbon_quota: float = None, --> carbon limitation (example : 99999999999999999999) @@ -87,39 +90,55 @@ step : int = 1, --> the number of step of budget where each time the candidates are registered time_index: int = 24*7, --> duration over which we optimize (in hour) time_interval: float = 1, --> the number of step of time where we evaluate the loss function (by default each 1 hour) - penalisation : float = None, --> penalisation cost + penalisation : float = None, --> penalisation cost carbon_cost : float = None, --> penalisation cost of carbon_quota plot : str = "default", --> write "None" for no plot average_wide : int = 0 --> average_wide of moving average plot parameter - - Output: + + Output: ------- List of results in each step of budget (here step = 20, that means only one step cause budget is set to 20) --> each result is a dict of "loss", "coef", "production", "unsatisfied demand", "carbon production" and "elapsed_time" """ -print(mix.optimizeMix(99999999999999999999,optimizer = opt_CMA, step = 20, penalisation = 100, carbon_cost = 0, time_index = 168, plot = "None"),"num_worker <------------ 1") -print(mix.optimizeMix(99999999999999999999,optimizer = opt_CMA_30, step = 20, penalisation = 100, carbon_cost = 0, time_index = 168, plot = "None"),"num_worker <------------ 30") -### Get all parameters used by the mix +print( + mix.optimizeMix(1e10, optimizer=opt_CMA, step=20, penalisation=100, carbon_cost=0, time_index=168, plot="None"), + "num_worker <------------ 1", +) +print( + mix.optimizeMix(1e10, optimizer=opt_CMA_30, step=20, penalisation=100, carbon_cost=0, time_index=168, plot="None"), + "num_worker <------------ 30", +) +### Get all parameters used by the mix print(mix.get_params()) -""" +""" (8) Evaluation of results by budget for each selected optimizer Default parameters : -------------------- mix, --> the mix to evaluate sequence, --> each budget to evaluate - max_budgets, - optimizer_list: List['str'], + max_budgets, + optimizer_list: List['str'], indicator_list: List['str'], --> indicators are ["loss","elapsed_time","production","unsatisfied demand","carbon production"] - num_worker: int = 1, + num_worker: int = 1, bind: str = None, --> path to dataset - time_index: int = 24, + time_index: int = 24, carbonProdLimit: float = 39500000000, --> equal to carbon_quota - time_interval : float = 1, - average_wide : int = 0, + time_interval : float = 1, + average_wide : int = 0, penalisation : float = 1000000000000, --> equal to penalisation cost carbon_cost: float = 0 """ -eva=EvaluationBudget() -eva.evaluate(mix,10,1000,optimizer_list = ["OnePlusOne","DE","CMA","PSO","NGOpt"], indicator_list = ["loss","elapsed_time","production","unsatisfied demand","carbon production"],carbonProdLimit = 9999999999999, time_index = 12, penalisation = 100, carbon_cost = 10) +eva = EvaluationBudget() +eva.evaluate( + mix, + 10, + 100, + optimizer_list=["OnePlusOne", "DE", "CMA", "PSO", "NGOpt"], + indicator_list=["loss", "elapsed_time", "production", "unsatisfied demand", "carbon production"], + carbonProdLimit=1e10, + time_index=24, + penalisation=100, + carbon_cost=10, +) diff --git a/mixsimulator/test_mas.py b/mixsimulator/test_mas.py new file mode 100644 index 0000000..0d5422b --- /dev/null +++ b/mixsimulator/test_mas.py @@ -0,0 +1,191 @@ +import random +import threading +import time +from typing import Any, Dict, List + +import numpy # type: ignore + +import mixsimulator.nevergradBased.Optimizer as opt +from mixsimulator import ElectricityMix +from mixsimulator.agents.Moderator import StoppableThread +from mixsimulator.Evaluation import EvaluationBudget + +""" + (0) Check the thread running in background +""" + + +def generate_random_scenario(centrals: List, time_index: int) -> Dict: + scenario = {} + for central in centrals: + tmp: Dict[Any, Any] = {"down": [], "up": []} + default_proba = random.uniform(0, 0.2) + + for i in range(time_index): + tmp[numpy.random.choice(["up", "down"], p=[1 - default_proba, default_proba])].append(i) + + up = [] + down = [] + for i in range(time_index): + if i in tmp["down"] and (i - 1 not in tmp["down"] or i == 0): + down.append(i) + if i - 1 in tmp["down"] and i in tmp["up"]: + up.append(i) + tmp["up"] = up + tmp["down"] = down + + scenario.update({central: tmp.copy()}) + + event_stack: Dict[Any, Any] = {} + for i in range(time_index): + for central in scenario.keys(): + if i in scenario[central]["down"]: + try: + event_stack[i].append(central._notify_is_down) + except: + event_stack.update({i: [central._notify_is_down]}) + elif i in scenario[central]["up"]: + try: + event_stack[i].append(central._notify_is_up) + except: + event_stack.update({i: [central._notify_is_up]}) + + # print(numpy.arange(0, 2)) + print("scenario: ", event_stack) + return event_stack + + +def check_thread_running(): + list_ = [] + while True: + tmp = threading.enumerate().copy() + if tmp != list_: + list_ = tmp + for thread in list_: + if thread.is_alive(): + print("THREAD: " + thread.name) + time.sleep(10) + + +thread_checker = StoppableThread(target=check_thread_running, name="thread_checker") +thread_checker.start() + +""" +(1) Configure nevergrad optimizers + + Default Parameters : + ---------- + opt = [ng.optimizers.OnePlusOne], + budget: List[int] = [20], + num_worker: int = 1, + instrum = ng.p.Array(shape=(2,)) +""" +opt_OPO = opt.Optimizer(opt=["OnePlusOne"], budget=[20], num_worker=1) +opt_OPO_20 = opt.Optimizer(opt=["OnePlusOne"], budget=[20], num_worker=20) + +""" +(2) Init MixSimulator instance : + Case one [Default] : "classic" method (see test_classic.py for more use case) + Case two : "MAS" or Multi Agent System method + + Default parameters : + ------------------------ + method : string = "classic", --> method explain above + carbon_cost : float = 0 --> cost of the CO2 production + penalisation_cost: float = 1e7 --> penalisation cost when production is more or less than the demand #NEED VERIFICATION +""" + +classic_mix = ElectricityMix().mix(method="classic", carbon_cost=0, penalisation_cost=100) +mas_mix = ElectricityMix().mix(method="MAS", carbon_cost=0, penalisation_cost=100) + +""" +(7) ONE SHOT optimization by calling the moderator of the MAS platform + + Default parameters : + --------------------------- + carbon_quota: float = None, --> carbon limitation (example : 99999999999999999999) + demand: Demand = None, --> you can use it to fix a constant demand for each time step (eliminates (5)) + lost: float = None, --> you can use it to fix a constant lost for each time step + optimizer: Optimizer = None, --> nevergrad optimizer initiate in (1) + step : int = 1, --> the number of step of budget where each time the candidates are registered + time_index: int = 24*7, --> duration over which we optimize (in hour) + time_interval: float = 1, --> the number of step of time where we evaluate the loss function (by default each 1 hour) + penalisation : float = None, --> penalisation cost + carbon_cost : float = None, --> penalisation cost of carbon_quota + plot : str = "default", --> write "None" for no plot + average_wide : int = 0 --> average_wide of moving average plot parameter + + Output: + ------- + List of results in each step of budget (here step = 20, that means only one step cause budget is set to 20) + --> each result is a dict of "loss", "coef", "production", "unsatisfied demand", "carbon production" and "elapsed_time" + +""" +print( + mas_mix.get_moderator().optimizeMix( + 1e10, optimizer=opt_OPO, step=20, penalisation=100, carbon_cost=0, time_index=168, plot="default" + ) +) + + +""" +(8) Evaluation of results by budget for each selected optimizer + Default parameters : + -------------------- + mix (or moderator), --> the mix or moderator to evaluate + sequence, --> each budget to evaluate + max_budgets, + optimizer_list: List['str'], + indicator_list: List['str'], --> indicators are ["loss","elapsed_time","production","unsatisfied demand","carbon production"] + num_worker: int = 1, + bind: str = None, --> path to dataset + time_index: int = 24, + carbonProdLimit: float = 39500000000, --> equal to carbon_quota + time_interval : float = 1, + average_wide : int = 0, + penalisation : float = 1000000000000, --> equal to penalisation cost + carbon_cost: float = 0 +""" +eva = EvaluationBudget() +eva.evaluate( + mas_mix.get_moderator(), + 10, + 100, + optimizer_list=["OnePlusOne", "DE", "CMA", "PSO", "NGOpt"], + indicator_list=["loss", "elapsed_time", "production", "unsatisfied demand", "carbon production"], + carbonProdLimit=1e10, + time_index=24, + penalisation=100, + carbon_cost=10, +) + +""" +(9) Simulating the mas platform (Manually) + 1 - First, set params by using set_params method + 2 - Run the run_optimization method to initiate the simulation + 3 - Add events +""" +mas_mix.get_moderator().set_params( + 1e10, optimizer=opt_OPO_20, step=20, penalisation=100, carbon_cost=0, time_index=4, plot="save" +) +mas_mix.get_moderator().run_optimization() + +# centrale1 = mas_mix.get_moderator().get_observable()[0] +# centrale2 = mas_mix.get_moderator().get_observable()[1] + +# centrale1._notify_is_down(6) +# centrale1._notify_is_up(12) +centrals = mas_mix.get_moderator().get_observable() +scenario = generate_random_scenario(centrals, 4) +for t in scenario.keys(): + for event in scenario[t]: + event(t) + +while True: + if len(threading.enumerate()) == 2: + thread_checker.stop() + break +print("SIMULATION DONE") + +print("FINAL RESULT: ", mas_mix.get_moderator().get_results()) +mas_mix.get_moderator().plotResults(mas_mix.get_moderator().get_results()) diff --git a/nevergradBased/Optimizer.py b/nevergradBased/Optimizer.py deleted file mode 100644 index 146fbf3..0000000 --- a/nevergradBased/Optimizer.py +++ /dev/null @@ -1,124 +0,0 @@ -import nevergrad as ng -import numpy as np # type: ignore -import math -import time -from typing import Type, List - -class Optimizer(): - """ - adaptation of nevergrad optimizers to the project + auto-parametrization - - list of available: self.getOptimizerList() - - """ - def __init__(self, opt = [ng.optimizers.OnePlusOne], budget: List[int] = [100], num_worker: int = 1, instrum = ng.p.Array(shape=(2,))): - self.set_budget(budget) - self.set_parametrization(instrum) - self.set_num_worker(num_worker) - - ### available optimizers - self.__available_optimizers = {Type[str]:Type[object]} - nevergrad_optimizers = list(ng.optimizers.registry.keys()) - for ng_opt in nevergrad_optimizers: - self.__available_optimizers.update({ng_opt:ng.optimizers.registry[ng_opt]}) - - convert_opt = [] - for opt_ng in opt: - if opt_ng in list(self.__available_optimizers.values()): - convert_opt.append(opt_ng) - continue - try: - convert_opt.append(self.__available_optimizers[str(opt_ng)]) - except: - raise KeyError(str(opt_ng) + " not available") - self.__optimizers = convert_opt - - def get_available_optimizers(self): - tmp = Optimizer() - return tmp.__available_optimizers.keys() - - def get_unavailable_optimizers(self): - tmp = Optimizer() - result = [] - tmp = tmp.__available_optimizers.keys() - list = sorted(ng.optimizers.registry.keys()) - for opt in list: - if opt not in tmp: - result.append(opt) - return result - - def get_optimizers(self): - return self.__optimizers - - def set_parametrization(self, instrum): - self.__parametrization = instrum - - def get_parametrization(self): - return self.__parametrization - - def set_budget(self,budget: List[int] = [100]): - #setting budget - #possible to do? auto calculate a optimal budget depending of the size of data - self.__budget = budget - - def get_budget(self): - return self.__budget - - def set_num_worker(self,n_work): - #setting budget - #possible to do? auto calculate a optimal budget depending of the size of data - self.__num_worker = n_work - - def get_num_worker(self): - return self.__num_worker - - def get_params(self) -> dict: - return {"optimizer(s)" : self.__optimizers, "budget(s)" : self.get_budget(), - "parametrization" : self.get_parametrization(), "num_worker" : self.get_num_worker()} - - def optimize(self, mix = None , func_to_optimize = None, constraints = None, step : int = 1): - - #setting budgets - budgets = self.get_budget() - - if len(budgets) != len(self.__optimizers) : - raise IndexError ("\n The length of budgets and the length of optimizers list should be the same.\n") - - total_budget = 0 - for algo_budget in budgets: - total_budget += algo_budget - result = [] - start_time = time.time() - - #optimization under constraints - chaining_algo = ng.optimizers.Chaining(self.__optimizers, budgets[:-1]) - optimizer = chaining_algo(parametrization=self.get_parametrization(), budget = budgets[-1], num_workers=self.get_num_worker()) - - #let's minimize - for tmp_budget in range(0, total_budget): - x = optimizer.ask() - loss = func_to_optimize(*x.args, constraints["time_interval"]) - optimizer.tell(x, loss) - if (tmp_budget+1)%step == 0: - result_per_budget = {} - recommendation = optimizer.provide_recommendation() - result_per_budget.update({"loss": func_to_optimize(recommendation.value, constraints["time_interval"])}) - result_per_budget.update({"coef": recommendation.value}) - - #Readjustment - if mix is not None : - usage_coef = mix.arrange_coef_as_array_of_array(recommendation.value) - weighted_coef = mix.get_weighted_coef(usage_coef, time_interval=constraints["time_interval"]) - production = 0 - u_demand = 0 - for t in range(0, len(weighted_coef)): - production += mix.get_production_cost_at_t(weighted_coef[t], t, constraints["time_interval"]) - u_demand += mix.get_unsatisfied_demand_at_t(weighted_coef[t], t, constraints["time_interval"]) - result_per_budget.update({"production": production}) - result_per_budget.update({"unsatisfied demand": u_demand}) - result_per_budget.update({"carbon production": mix.get_carbon_production_at_t(weighted_coef[t], constraints["time_interval"])}) - result_per_budget.update({"elapsed_time": time.time() - start_time}) - result.append(result_per_budget) - - #TODO --> Check Constraints - return result diff --git a/pyproject.toml b/pyproject.toml new file mode 100644 index 0000000..ba17c85 --- /dev/null +++ b/pyproject.toml @@ -0,0 +1,26 @@ + +# Example configuration for Black. + +# NOTE: you have to use single-quoted strings in TOML for regular expressions. +# It's the equivalent of r-strings in Python. Multiline strings are treated as +# verbose regular expressions by Black. Use [ ] to denote a significant space +# character. + +[tool.black] +line-length = 119 +target-version = ['py36', 'py37', 'py38'] +include = '\.pyi?$' +exclude = ''' +/( + \.eggs + | \.git + | \.hg + | \.mypy_cache + | \.tox + | \.venv + | _build + | buck-out + | build + | dist +)/ +''' diff --git a/requirements_dev.txt b/requirements_dev.txt new file mode 100644 index 0000000..0d7c87e --- /dev/null +++ b/requirements_dev.txt @@ -0,0 +1,20 @@ +black==21.12b0 +flake8==4.0.1 +iglovikov-helper-functions==0.0.53 +isort==5.10.1 +mypy==0.902 +pep8-naming==0.11.1 +pip==21.3.1 +pre_commit==2.17.0 +pylint==2.12.2 +pytest==6.2.4 +pytest-asyncio==0.15.1 +pytest-clarity==0.3.0a0 +pytest-cov==2.11.1 +pytest-mock==3.6.1 +pytest-parallel==0.1.0 +pytest-sugar==0.9.4 +pytest-xdist==2.2.1 +types-PyYAML==5.4.3 +types-requests==2.25.0 + diff --git a/setup.cfg b/setup.cfg new file mode 100644 index 0000000..e34c821 --- /dev/null +++ b/setup.cfg @@ -0,0 +1,15 @@ +[flake8] +max-line-length = 120 +exclude =.git,__pycache__,docs/source/conf.py,build,dist,tests +ignore = I101,D100,D101,D102,D103,D104,D105,D107,D401,E203,I900,N802,N806,N812,W503,S311,S605,S607 + +[mypy] +ignore_missing_imports = True +disallow_untyped_defs = True +check_untyped_defs = True +warn_redundant_casts = True +no_implicit_optional = True +strict_optional = True + +[mypy-tests.*] +ignore_errors = True diff --git a/setup.py b/setup.py new file mode 100644 index 0000000..f5cc542 --- /dev/null +++ b/setup.py @@ -0,0 +1,35 @@ +import setuptools # type: ignore + +with open("README.md") as fh: + long_description = fh.read() + +setuptools.setup( + name="mixsimulator", + version="0.4.5", + author="RASOANAIVO Andry, ANDRIAMALALA Rahamefy Solofohanitra, ANDRIAMIZAKASON Toky Axel", + author_email="tokyandriaxel@gmail.com", + description="Python application with nevergrad optimization model for calculating and simulating the least cost of an energy Mix under constraints.", + long_description=long_description, + long_description_content_type="text/markdown", + url="https://github.com/Foloso/MixSimulator", + packages=setuptools.find_packages(include=["mixsimulator", "mixsimulator.*"]), + include_package_data=True, + package_data={ + "mixsimulator": [ + "Experiments/Scenario_type.py", + "data/RIToamasina/dataset_RI_Toamasina.csv", + "data/RIToamasina/dataset_RI_Toamasina_v2.csv", + "data/RIToamasina/dataset_RI_Toamasina_variation_template.csv", + "data/RIToamasina/DIR-TOAMASINA_concat.csv", + "LICENSE", + "params_files/exchange_code.json", + "params_files/settings.json", + ] + }, + classifiers=[ + "Programming Language :: Python :: 3", + "License :: OSI Approved :: MIT License", + "Operating System :: OS Independent", + ], + python_requires=">=3.6", +)