random decider for debugging

v2/architecture.py

@@ -5,7 +5,7 @@ import time
 import matplotlib.pyplot as plt
 
 # it's really just a network pricing model
-from decider import Decider, Decision
+from decider import Decision, Decider, RandomDecider
 from supplier import Supplier, precalculate_supplier
 from data_processing import read_datasets, add_production_field, interpolate_and_join, SolarParameters
 from bess import BatteryModel
@@ -27,8 +27,8 @@ def add_dummy_predictions(all_data_with_predictions):
     all_data_with_predictions['Consumption_prediction'] = all_data_with_predictions['Consumption'].shift(periods=cons_shift)
     all_data_with_predictions['Production_prediction'] = all_data_with_predictions['Production'].shift(periods=prod_shift)
     # we predict zero before we have data, no big deal:
-    all_data_with_predictions['Consumption_prediction'][:cons_shift] = 0
-    all_data_with_predictions['Production_prediction'][:prod_shift] = 0
+    all_data_with_predictions.loc[all_data_with_predictions.index[:cons_shift], 'Consumption_prediction'] = 0
+    all_data_with_predictions.loc[all_data_with_predictions.index[:prod_shift], 'Production_prediction'] = 0
 
 
 # even mock-er class than usual.
@@ -185,7 +185,7 @@ def simulator(battery_model, prod_cons, decider):
     '''
     consumption_charge_series = consumption_fees_np * consumption_from_network_pandas_series.to_numpy()
     step_in_hour = consumption_from_network_pandas_series.index.freq.n / 60 # [hour], the length of a time step.
-    fifteen_minute_demands_in_kwh = consumption_from_network_pandas_series.resample('15T').sum() * step_in_hour
+    fifteen_minute_demands_in_kwh = consumption_from_network_pandas_series.resample('15min').sum() * step_in_hour
     fifteen_minute_surdemands_in_kwh = (fifteen_minute_demands_in_kwh - decider.precalculated_supplier.peak_demand).clip(lower=0)
     demand_charges = fifteen_minute_surdemands_in_kwh * decider.precalculated_supplier.surcharge_per_kwh
     total_network_fee = consumption_charge_series.sum() + demand_charges.sum()
@@ -210,23 +210,62 @@ def simulator(battery_model, prod_cons, decider):
 
 
 
-def optimizer(battery_model, all_data_with_predictions, precalculated_supplier):
+def optimizer(decider_class, battery_model, all_data_with_predictions, precalculated_supplier):
+    number_of_generations = 3
+    population_size = 100
+    collected_loss_values = []
     def objective_function(params):
         print("Simulating with parameters", params)
-        # TODO params completely ignored right now.
-        decider = Decider(params, precalculated_supplier)
+        decider = decider_class(params, precalculated_supplier)
         t = time.perf_counter()
         results, total_network_fee = simulator(battery_model, all_data_with_predictions, decider)
+        collected_loss_values.append((params, total_network_fee))
         return total_network_fee
 
     def clipper_function(params):
-        return Decider.clip_params(params)
+        return decider_class.clip_params(params)
+
+    init_mean, init_scale = decider_class.initial_params()
+    best_params = evolution_strategies_optimizer(objective_function, clipper_function,
+                                                 init_mean=init_mean, init_scale=init_scale,
+                                                 number_of_generations=number_of_generations,
+                                                 population_size=population_size)
+    return best_params, collected_loss_values
+
+
+def visualize_collected_loss_values(collected_loss_values):
+    all_params = []
+    losses = []
+    for row in collected_loss_values:
+        params, loss = row
+        all_params.append(params)
+        losses.append(loss)
+
+    all_params = np.array(all_params)
+    losses = np.array(losses)
+
+    # losses -= losses.min() ; losses /= losses.max()
+
+    plt.scatter(all_params[:, 0], all_params[:, 1], c=range(len(all_params)))
+    plt.show()
+
+
+    from mpl_toolkits.mplot3d import Axes3D
+
+    fig = plt.figure()
+
+    # Add a 3D subplot
+    ax = fig.add_subplot(111, projection='3d')
+
+    # Scatter plot
+    ax.scatter(all_params[:, 0], all_params[:, 1], losses)
+    plt.show()
 
-    init_mean, init_scale = Decider.initial_params()
-    best_params = evolution_strategies_optimizer(objective_function, clipper_function, init_mean=init_mean, init_scale=init_scale)
 
 
 def main():
+    np.random.seed(1)
+
     supplier = Supplier(price=100) # Ft/kWh
     # nine-to-five increased price.
     supplier.set_price_for_daily_interval(9, 17, 150)
@@ -242,11 +281,17 @@ def main():
     add_production_field(met_2021_data, solar_parameters)
     all_data = interpolate_and_join(met_2021_data, cons_2021_data)
 
-    print("Working with", solar_parameters.solar_cell_num, "solar cells, that's a maximum production of", all_data['Production'].max(), "kW.")
-
     time_interval_min = all_data.index.freq.n
     time_interval_h = time_interval_min / 60
 
+    # for faster testing:
+    DATASET_TRUNCATED_SIZE = 10000
+    if DATASET_TRUNCATED_SIZE is not None:
+        print("Truncating dataset to", DATASET_TRUNCATED_SIZE, "datapoints, that is", DATASET_TRUNCATED_SIZE * time_interval_h / 24, "days")
+        all_data = all_data.iloc[:DATASET_TRUNCATED_SIZE]
+
+    print("Working with", solar_parameters.solar_cell_num, "solar cells, that's a maximum production of", all_data['Production'].max(), "kW.")
+
     all_data_with_predictions = all_data.copy()
     add_dummy_predictions(all_data_with_predictions)
 
@@ -258,11 +303,13 @@ def main():
 
     battery_model = BatteryModel(capacity_Ah=600, time_interval_h=time_interval_h)
 
+    decider_class = RandomDecider
+
     # TODO this is super unfortunate:
     # Consumption_fees travels via all_data_with_predictions,
     # peak_demand and surcharge_per_kwh travels via precalculated_supplier of decider.
-    decider_init_mean, decider_init_scale = Decider.initial_params()
-    decider = Decider(decider_init_mean, precalculated_supplier)
+    decider_init_mean, decider_init_scale = decider_class.initial_params()
+    decider = decider_class(decider_init_mean, precalculated_supplier)
 
     t = time.perf_counter()
     results, total_network_fee = simulator(battery_model, all_data_with_predictions, decider)
@@ -273,7 +320,11 @@ def main():
         plt.title('soc_series')
         plt.show()
 
-    optimizer(battery_model, all_data_with_predictions, precalculated_supplier)
+    best_params, collected_loss_values = optimizer(decider_class, battery_model, all_data_with_predictions, precalculated_supplier)
+    visualize_collected_loss_values(collected_loss_values)
+
+    decider = decider_class(best_params, precalculated_supplier)
+    results, total_network_fee = simulator(battery_model, all_data_with_predictions, decider)
 
 if __name__ == '__main__':
     main()

v2/decider.py

@@ -25,12 +25,50 @@ class Decision(IntEnum):
     NETWORK_CHARGE = 2
 
 
+class RandomDecider:
+    def __init__(self, params, precalculated_supplier):
+        self.input_window_size = STEPS_PER_HOUR * 24 # day long window.
+        self.precalculated_supplier = precalculated_supplier
+        assert params.shape == (2, )
+        # param_1 is prob of choosing PASSIVE
+        # param_2 is prob of choosing NETWORK_CHARGE
+        self.passive_prob, self.network_prob = params
+
+    def decide(self, prod_pred, cons_pred, fees, battery_model):
+        r = np.random.rand()
+        if r < self.passive_prob:
+            return Decision.PASSIVE
+        elif r < self.passive_prob + self.network_prob:
+            return Decision.NETWORK_CHARGE
+        else:
+            return Decision.DISCHARGE
+
+    @staticmethod
+    def clip_params(params):
+        assert params.shape == (2, )
+        p1, p2 = params
+        p1 = max((0, p1))
+        p2 = max((0, p2))
+        s = p1 + p2
+        if s > 1:
+            p1 /= s
+            p2 /= s
+        return np.array([p1, p2])
+
+
+    @staticmethod
+    def initial_params():
+        init_mean = np.array([1/3, 1/3])
+        init_scale = np.array([1.0, 1.0])
+        return init_mean, init_scale
+
+
+
 # mock class as usual
 # output_window_size is not yet used, always decides one timestep.
 class Decider:
     def __init__(self, params, precalculated_supplier):
         self.input_window_size = STEPS_PER_HOUR * 24 # day long window.
-        self.random_seed = 0
         self.precalculated_supplier = precalculated_supplier
         assert params.shape == (2, )
         # param_1 is how many minutes do we look ahead to decide if there's

v2/evolution_strategies.py

@@ -11,10 +11,11 @@ def clip_params(population, clipper_function):
         population[i] = clipper_function(individual)
 
 
-def evolution_strategies_optimizer(objective_function, clipper_function, init_mean, init_scale):
+def evolution_strategies_optimizer(objective_function, clipper_function,
+                                   init_mean, init_scale,
+                                   number_of_generations,
+                                   population_size):
     # Initialize parameters
-    population_size = 100
-    number_of_generations = 30
     mutation_scale = 0.05
     selection_ratio = 0.5
     selected_size = int(population_size * selection_ratio)
@@ -44,19 +45,6 @@ def evolution_strategies_optimizer(objective_function, clipper_function, init_me
         best_index = np.argmin(fitness)
         best_solution = population[best_index]
         print(f"Generation {generation + 1}: Best Fitness = {best_fitness}", f"Best solution so far: {best_solution}")
-        if DO_VIS:
-            plt.scatter(population[:, 0], population[:, 1])
-            plt.xlim(-5, 5)
-            plt.ylim(-5, 5)
-            x = np.linspace(-5, 5, 100)
-            y = np.linspace(-5, 5, 100)
-            XY = np.stack(np.meshgrid(x, y), axis=-1)
-            print(XY.shape)
-
-            # Defining the function (x - 3)**2 + (y + 2)**2
-            Z = toy_objective_function(XY)
-            contour = plt.contour(XY[..., 0], XY[..., 1], Z, levels=50)
-            plt.show()
 
 
     # Best solution

v2/supplier.py

@@ -71,7 +71,7 @@ class Supplier:
         assert 15 % demand_series.index.freq.n == 0
         time_steps_per_demand_charge_evaluation = 15 // demand_series.index.freq.n
         # fifteen_minute_peaks [kW] tells the maximum demand in a 15 minutes timeframe:
-        fifteen_minute_demands_in_kwh = demand_series.resample('15T').sum() * step_in_hour
+        fifteen_minute_demands_in_kwh = demand_series.resample('15min').sum() * step_in_hour
         demand_charges = pd.Series([self.demand_charge(demand_in_kwh) for demand_in_kwh in fifteen_minute_demands_in_kwh], index=fifteen_minute_demands_in_kwh.index)
         total_demand_charge = sum(demand_charges)
         total_charge = consumption_charge + total_demand_charge