Merge pull request #142 from davidusb-geek/davidusb-geek/fix/constrai…

…nts_conflict

Fix - Fixed conflict on constraints

CHANGELOG.md

@@ -1,5 +1,12 @@
 # Changelog
 
+## [0.6.2] - 2024-01-04
+### Improvement
+- Added option to pass additional weight for battery usage
+- Improved coverage
+### Fix
+- Updated optimization constraints to solve conflict for `set_def_constant` and `treat_def_as_semi_cont` cases
+
 ## [0.6.1] - 2023-12-18
 ### Fix
 - Patching EMHASS for Python 3.11. New explicit dependecy h5py==3.10.0

Dockerfile

@@ -1,5 +1,6 @@
 FROM python:3.11-slim-buster
-#FROM ghcr.io/home-assistant/amd64-base-debian:bookworm # Uncomment to test add-on
+# FROM ghcr.io/home-assistant/amd64-base-debian:bookworm # Uncomment to test add-on
+# FROM ghcr.io/home-assistant/armhf-base-debian:bookworm
 
 # switch working directory
 WORKDIR /app
@@ -13,22 +14,25 @@ COPY README.md README.md
 # Setup
 RUN apt-get update \
     && apt-get install -y --no-install-recommends \
-        # libc-bin \ # Uncomment to test add-on
-        # libffi-dev \ # Uncomment to test add-on
-        # python3 \ # Uncomment to test add-on
-        # python3-pip \ # Uncomment to test add-on
-        # python3-dev \ # Uncomment to test add-on
-        # git \ # Uncomment to test add-on
-        # build-essential \ # Uncomment to test add-on
+        # libffi-dev \
+        # python3 \
+        # python3-pip \
+        # python3-dev \
+        # git \
+        # build-essential \
         gcc \
         coinor-cbc \
         coinor-libcbc-dev \
         libhdf5-dev \
         libhdf5-serial-dev \
         netcdf-bin \
         libnetcdf-dev \
+        # pkg-config \
+        # gfortran \
+        # libatlas-base-dev \
     && ln -s /usr/include/hdf5/serial /usr/include/hdf5/include \
     && export HDF5_DIR=/usr/include/hdf5 \
+    # && pip3 install --extra-index-url=https://www.piwheels.org/simple --no-cache-dir --break-system-packages -U setuptools wheel \
     && pip3 install --no-cache-dir --break-system-packages -r requirements_webserver.txt \
     && apt-get purge -y --auto-remove \
         gcc \

config_emhass.yaml

@@ -52,6 +52,8 @@ optim_conf:
   - set_battery_dynamic: False # add a constraint to limit the dynamic of the battery power in power per time unit
   - battery_dynamic_max: 0.9 # maximum dynamic positive power variation in percentage of battery maximum power
   - battery_dynamic_min: -0.9 # minimum dynamic negative power variation in percentage of battery maximum power
+  - weight_battery_discharge: 1.0 # weight applied in cost function to battery usage for discharge
+  - weight_battery_charge: 1.0 # weight applied in cost function to battery usage for charge
 
 plant_conf:
   - P_grid_max: 9000 # The maximum power that can be supplied by the utility grid in Watts

docs/conf.py

@@ -22,7 +22,7 @@
 author = 'David HERNANDEZ'
 
 # The full version, including alpha/beta/rc tags
-release = '0.6.1'
+release = '0.6.2'
 
 # -- General configuration ---------------------------------------------------
 

docs/config.md

@@ -79,6 +79,8 @@ The following parameters and definitions are only needed if load_cost_forecast_m
 - set_battery_dynamic: Set a power dynamic limiting condition to the battery power. This is an additional constraint on the battery dynamic in power per unit of time, which allows you to set a percentage of the battery nominal full power as the maximum power allowed for (dis)charge.
 - battery_dynamic_max: The maximum positive (for discharge) battery power dynamic. This is the allowed power variation (in percentage) of battery maximum power per unit of time.
 - battery_dynamic_min: The maximum negative (for charge) battery power dynamic. This is the allowed power variation (in percentage) of battery maximum power per unit of time.
+- weight_battery_discharge: An additional weight applied in cost function to battery usage for discharge.
+- weight_battery_charge: An additional weight applied in cost function to battery usage for charge.
 
 ## System configuration parameters
 

options.json

@@ -13,6 +13,8 @@
   "set_battery_dynamic": false,
   "battery_dynamic_max": 0.9,
   "battery_dynamic_min": -0.9,
+  "weight_battery_discharge": 1.0,
+  "weight_battery_charge": 1.0,
   "sensor_power_photovoltaics": "sensor.power_photovoltaics",                                                                                                                                                    
   "sensor_power_load_no_var_loads": "sensor.power_load_no_var_loads",                                                                                                                                            
   "number_of_deferrable_loads": 3,                                                                                                                                                                               

requirements.txt

@@ -1,7 +1,7 @@
 wheel
-numpy==1.26.0
-pandas==2.0.3
-scipy==1.11.3
+numpy<=1.26.0
+pandas<=2.0.3
+scipy<=1.11.3
 pvlib>=0.10.2
 protobuf>=3.0.0
 pytz>=2021.1
@@ -10,5 +10,5 @@ beautifulsoup4>=4.9.3
 h5py==3.10.0
 pulp>=2.4
 pyyaml>=5.4.1
-tables==3.9.1
+tables<=3.9.1
 skforecast==0.11.0

scripts/optim_results_analysis.py

@@ -0,0 +1,119 @@
+# -*- coding: utf-8 -*-
+import pickle
+import numpy as np
+import pandas as pd
+import pathlib
+import plotly.express as px
+import plotly.subplots as sp
+import plotly.io as pio
+pio.renderers.default = 'browser'
+pd.options.plotting.backend = "plotly"
+
+from emhass.retrieve_hass import retrieve_hass
+from emhass.optimization import optimization
+from emhass.forecast import forecast
+from emhass.utils import get_root, get_yaml_parse, get_days_list, get_logger
+
+# the root folder
+root = str(get_root(__file__, num_parent=2))
+# create logger
+logger, ch = get_logger(__name__, root, save_to_file=False)
+
+def get_forecast_optim_objects(retrieve_hass_conf, optim_conf, plant_conf,
+                               params, get_data_from_file):
+    fcst = forecast(retrieve_hass_conf, optim_conf, plant_conf,
+                    params, root, logger, get_data_from_file=get_data_from_file)
+    df_weather = fcst.get_weather_forecast(method='csv')
+    P_PV_forecast = fcst.get_power_from_weather(df_weather)
+    P_load_forecast = fcst.get_load_forecast(method=optim_conf['load_forecast_method'])
+    df_input_data_dayahead = pd.concat([P_PV_forecast, P_load_forecast], axis=1)
+    df_input_data_dayahead.columns = ['P_PV_forecast', 'P_load_forecast']
+    opt = optimization(retrieve_hass_conf, optim_conf, plant_conf, 
+                       fcst.var_load_cost, fcst.var_prod_price,  
+                       'profit', root, logger)
+    return fcst, P_PV_forecast, P_load_forecast, df_input_data_dayahead, opt
+
+if __name__ == '__main__':
+    show_figures = False
+    save_figures = False
+    get_data_from_file = True
+    params = None
+    retrieve_hass_conf, optim_conf, plant_conf = get_yaml_parse(pathlib.Path(root+'/config_emhass.yaml'), use_secrets=False)
+    retrieve_hass_conf, optim_conf, plant_conf = \
+        retrieve_hass_conf, optim_conf, plant_conf
+    rh = retrieve_hass(retrieve_hass_conf['hass_url'], retrieve_hass_conf['long_lived_token'], 
+                        retrieve_hass_conf['freq'], retrieve_hass_conf['time_zone'],
+                        params, root, logger)
+    if get_data_from_file:
+        with open(pathlib.Path(root+'/data/test_df_final.pkl'), 'rb') as inp:
+            rh.df_final, days_list, var_list = pickle.load(inp)
+    else:
+        days_list = get_days_list(retrieve_hass_conf['days_to_retrieve'])
+        var_list = [retrieve_hass_conf['var_load'], retrieve_hass_conf['var_PV']]
+        rh.get_data(days_list, var_list,
+                        minimal_response=False, significant_changes_only=False)
+    rh.prepare_data(retrieve_hass_conf['var_load'], load_negative = retrieve_hass_conf['load_negative'],
+                            set_zero_min = retrieve_hass_conf['set_zero_min'], 
+                            var_replace_zero = retrieve_hass_conf['var_replace_zero'], 
+                            var_interp = retrieve_hass_conf['var_interp'])
+    df_input_data = rh.df_final.copy()
+
+    fcst, P_PV_forecast, P_load_forecast, df_input_data_dayahead, opt = \
+        get_forecast_optim_objects(retrieve_hass_conf, optim_conf, plant_conf,
+                                   params, get_data_from_file)
+    df_input_data = fcst.get_load_cost_forecast(df_input_data)
+    df_input_data = fcst.get_prod_price_forecast(df_input_data)
+
+    template = 'presentation'
+
+    # Let's plot the input data
+    fig_inputs1 = df_input_data[['sensor.power_photovoltaics',
+                                'sensor.power_load_no_var_loads_positive']].plot()
+    fig_inputs1.layout.template = template
+    fig_inputs1.update_yaxes(title_text = "Powers (W)")
+    fig_inputs1.update_xaxes(title_text = "Time")
+    if show_figures:
+        fig_inputs1.show()
+    if save_figures:
+        fig_inputs1.write_image(root + "/docs/images/inputs_power.svg", 
+                                width=1080, height=0.8*1080)
+
+    fig_inputs_dah = df_input_data_dayahead.plot()
+    fig_inputs_dah.layout.template = template
+    fig_inputs_dah.update_yaxes(title_text = "Powers (W)")
+    fig_inputs_dah.update_xaxes(title_text = "Time")
+    if show_figures:
+        fig_inputs_dah.show()
+    if save_figures:
+        fig_inputs_dah.write_image(root + "/docs/images/inputs_dayahead.svg", 
+                                   width=1080, height=0.8*1080)
+
+    # And then perform a dayahead optimization
+    df_input_data_dayahead = fcst.get_load_cost_forecast(df_input_data_dayahead)
+    df_input_data_dayahead = fcst.get_prod_price_forecast(df_input_data_dayahead)
+    optim_conf['treat_def_as_semi_cont'] = [True, True]
+    optim_conf['set_def_constant'] = [True, True]
+    unit_load_cost = df_input_data[opt.var_load_cost].values
+    unit_prod_price = df_input_data[opt.var_prod_price].values
+    opt_res_dah = opt.perform_optimization(df_input_data_dayahead, P_PV_forecast.values.ravel(), 
+                                           P_load_forecast.values.ravel(), 
+                                           unit_load_cost, unit_prod_price,
+                                           debug = True)
+    # opt_res_dah = opt.perform_dayahead_forecast_optim(df_input_data_dayahead, P_PV_forecast, P_load_forecast)
+    opt_res_dah['P_PV'] = df_input_data_dayahead[['P_PV_forecast']]
+    fig_res_dah = opt_res_dah[['P_deferrable0', 'P_deferrable1', 'P_grid', 'P_PV',
+                               'P_def_start_0', 'P_def_start_1', 'P_def_bin2_0', 'P_def_bin2_1']].plot()
+    fig_res_dah.layout.template = template
+    fig_res_dah.update_yaxes(title_text = "Powers (W)")
+    fig_res_dah.update_xaxes(title_text = "Time")
+    # if show_figures:
+    fig_res_dah.show()
+    if save_figures:
+        fig_res_dah.write_image(root + "/docs/images/optim_results_PV_defLoads_dayaheadOptim.svg", 
+                                width=1080, height=0.8*1080)
+
+    print("System with: PV, two deferrable loads, dayahead optimization, profit >> total cost function sum: "+\
+        str(opt_res_dah['cost_profit'].sum()))
+
+    print(opt_res_dah)
+    opt_res_dah.to_html('opt_res_dah.html')

setup.py

@@ -19,7 +19,7 @@
 
 setup(
     name='emhass',  # Required
-    version='0.6.1',  # Required
+    version='0.6.2',  # Required
     description='An Energy Management System for Home Assistant',  # Optional
     long_description=long_description,  # Optional
     long_description_content_type='text/markdown',  # Optional (see note above)
@@ -40,17 +40,17 @@
     python_requires='>=3.9, <3.12',
     install_requires=[
         'wheel', 
-        'numpy==1.26',
-        'scipy==1.11.3',
-        'pandas==2.0.3',
+        'numpy<=1.26',
+        'scipy<=1.11.3',
+        'pandas<=2.0.3',
         'pvlib>=0.10.1',
         'protobuf>=3.0.0',
         'pytz>=2021.1',
         'requests>=2.25.1',
         'beautifulsoup4>=4.9.3',
         'pulp>=2.4',
         'pyyaml>=5.4.1',
-        'tables==3.9.1',
+        'tables<=3.9.1',
         'skforecast==0.11.0',
     ],  # Optional
     entry_points={  # Optional

src/emhass/optimization.py

@@ -77,7 +77,7 @@ def __init__(self, retrieve_hass_conf: dict, optim_conf: dict, plant_conf: dict,
         if 'lp_solver' in optim_conf.keys():
             self.lp_solver = optim_conf['lp_solver']
         else:
-            self.lp_solver = 'PULP_CBC_CMD'
+            self.lp_solver = 'default'
         if 'lp_solver_path' in optim_conf.keys():
             self.lp_solver_path = optim_conf['lp_solver_path']
         else:
@@ -88,7 +88,8 @@ def __init__(self, retrieve_hass_conf: dict, optim_conf: dict, plant_conf: dict,
     def perform_optimization(self, data_opt: pd.DataFrame, P_PV: np.array, P_load: np.array, 
                              unit_load_cost: np.array, unit_prod_price: np.array,
                              soc_init: Optional[float] = None, soc_final: Optional[float] = None,
-                             def_total_hours: Optional[list] = None) -> pd.DataFrame:
+                             def_total_hours: Optional[list] = None, 
+                             debug: Optional[bool] = False) -> pd.DataFrame:
         r"""
         Perform the actual optimization using linear programming (LP).
         
@@ -196,31 +197,31 @@ def perform_optimization(self, data_opt: pd.DataFrame, P_PV: np.array, P_load: n
         if self.costfun == 'profit':
             if self.optim_conf['set_total_pv_sell']:
                 objective = plp.lpSum(-0.001*self.timeStep*(unit_load_cost[i]*(P_load[i] + P_def_sum[i]) + \
-                                                            unit_prod_price[i]*P_grid_neg[i])
-                                      for i in set_I)
+                    unit_prod_price[i]*P_grid_neg[i]) for i in set_I)
             else:
                 objective = plp.lpSum(-0.001*self.timeStep*(unit_load_cost[i]*P_grid_pos[i] + \
-                                                            unit_prod_price[i]*P_grid_neg[i])
-                                      for i in set_I)
+                    unit_prod_price[i]*P_grid_neg[i]) for i in set_I)
         elif self.costfun == 'cost':
             if self.optim_conf['set_total_pv_sell']:
-                objective = plp.lpSum(-0.001*self.timeStep*unit_load_cost[i]*(P_load[i] + P_def_sum[i])
-                                      for i in set_I)
+                objective = plp.lpSum(-0.001*self.timeStep*unit_load_cost[i]*(P_load[i] + P_def_sum[i]) for i in set_I)
             else:
-                objective = plp.lpSum(-0.001*self.timeStep*unit_load_cost[i]*P_grid_pos[i]
-                                      for i in set_I)
+                objective = plp.lpSum(-0.001*self.timeStep*unit_load_cost[i]*P_grid_pos[i] for i in set_I)
         elif self.costfun == 'self-consumption':
             if type_self_conso == 'bigm':
                 bigm = 1e3
                 objective = plp.lpSum(-0.001*self.timeStep*(bigm*unit_load_cost[i]*P_grid_pos[i] + \
-                                                            unit_prod_price[i]*P_grid_neg[i])
-                                      for i in set_I)
+                    unit_prod_price[i]*P_grid_neg[i]) for i in set_I)
             elif type_self_conso == 'maxmin':
                 objective = plp.lpSum(0.001*self.timeStep*unit_load_cost[i]*SC[i] for i in set_I)
             else:
                 self.logger.error("Not a valida option for type_self_conso parameter")
         else:
             self.logger.error("The cost function specified type is not valid")
+        # Add more terms to the objective function in the case of battery use
+        if self.optim_conf['set_use_battery']:
+            objective = objective + plp.lpSum(-0.001*self.timeStep*(
+                self.optim_conf['weight_battery_discharge']*P_sto_pos[i] + \
+                    self.optim_conf['weight_battery_charge']*P_sto_neg[i]) for i in set_I)
         opt_model.setObjective(objective)
 
         ## Setting constraints
@@ -287,30 +288,31 @@ def perform_optimization(self, data_opt: pd.DataFrame, P_PV: np.array, P_load: n
                     for i in set_I})
             # Treat the number of starts for a deferrable load
             if self.optim_conf['set_def_constant'][k]:
-                constraints.update({"constraint_pdef{}_start1".format(k) : 
-                    plp.LpConstraint(
-                        e=P_def_start[k][0],
-                        sense=plp.LpConstraintEQ,
-                        rhs=0)
-                    })
-                constraints.update({"constraint_pdef{}_start2_{}".format(k, i) : 
-                    plp.LpConstraint(
-                        e=P_def_start[k][i] - P_def_bin2[k][i] + P_def_bin2[k][i-1],
-                        sense=plp.LpConstraintEQ,
-                        rhs=0)
-                    for i in set_I[1:]})
-                constraints.update({"constraint_pdef{}_start4_{}".format(k, i) : 
+
+                constraints.update({"constraint_pdef{}_start1_{}".format(k, i) : 
                     plp.LpConstraint(
                         e=P_deferrable[k][i] - P_def_bin2[k][i]*M,
                         sense=plp.LpConstraintLE,
                         rhs=0)
                     for i in set_I})
-                constraints.update({"constraint_pdef{}_start5_{}".format(k, i) : 
+                constraints.update({"constraint_pdef{}_start2_{}".format(k, i) : 
                     plp.LpConstraint(
-                        e=-P_deferrable[k][i] + M*(P_def_bin2[k][i]-1) + 1,
-                        sense=plp.LpConstraintLE,
+                        e=P_def_start[k][i] - P_def_bin2[k][i] + P_def_bin2[k][i-1],
+                        sense=plp.LpConstraintGE,
                         rhs=0)
-                    for i in set_I})
+                    for i in set_I[1:]})
+                constraints.update({"constraint_pdef{}_start3".format(k) :
+                plp.LpConstraint(
+                    e = plp.lpSum(P_def_start[k][i] for i in set_I),
+                    sense = plp.LpConstraintEQ,
+                    rhs = 1)
+                })
+                constraints.update({"constraint_pdef{}_start4".format(k) :
+                plp.LpConstraint(
+                    e = plp.lpSum(P_def_bin2[k][i] for i in set_I),
+                    sense = plp.LpConstraintEQ,
+                    rhs = self.optim_conf['def_total_hours'][k]/self.timeStep)
+                })
 
         # The battery constraints
         if self.optim_conf['set_use_battery']:
@@ -465,6 +467,13 @@ def perform_optimization(self, data_opt: pd.DataFrame, P_PV: np.array, P_load: n
         # Add the optimization status
         opt_tp["optim_status"] = self.optim_status
 
+        # Debug variables
+        if debug:
+            opt_tp["P_def_start_0"] = [P_def_start[0][i].varValue for i in set_I]
+            opt_tp["P_def_start_1"] = [P_def_start[1][i].varValue for i in set_I]
+            opt_tp["P_def_bin2_0"] = [P_def_bin2[0][i].varValue for i in set_I]
+            opt_tp["P_def_bin2_1"] = [P_def_bin2[1][i].varValue for i in set_I]
+
         return opt_tp
 
     def perform_perfect_forecast_optim(self, df_input_data: pd.DataFrame, days_list: pd.date_range) -> pd.DataFrame:

src/emhass/utils.py

@@ -324,6 +324,10 @@ def treat_runtimeparams(runtimeparams: str, params: str, retrieve_hass_conf: dic
         if 'solar_forecast_kwp' in runtimeparams.keys():
             retrieve_hass_conf['solar_forecast_kwp'] = runtimeparams['solar_forecast_kwp']
             optim_conf['weather_forecast_method'] = 'solar.forecast'
+        if 'weight_battery_discharge' in runtimeparams.keys():
+            optim_conf['weight_battery_discharge'] = runtimeparams['weight_battery_discharge']
+        if 'weight_battery_charge' in runtimeparams.keys():
+            optim_conf['weight_battery_charge'] = runtimeparams['weight_battery_charge']
         # Treat plant configuration parameters passed at runtime
         if 'SOCtarget' in runtimeparams.keys():
             plant_conf['SOCtarget'] = runtimeparams['SOCtarget']