diff --git a/example.ipynb b/example.ipynb index 879070a..68c7b5b 100644 --- a/example.ipynb +++ b/example.ipynb @@ -9,7 +9,7 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": 8, "metadata": {}, "outputs": [], "source": [ @@ -27,7 +27,7 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 9, "metadata": {}, "outputs": [], "source": [ @@ -64,16 +64,16 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ - "" + "" ] }, - "execution_count": 4, + "execution_count": 10, "metadata": {}, "output_type": "execute_result" }, @@ -117,14 +117,14 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "optimizer_flex = optimizer(solverpath_exe='C:\\glpk\\w64\\glpsol')\n", "step1_soc_daa,step1_cha_daa,step1_dis_daa, step1_profit_daa = optimizer_flex.step1_optimize_daa(n_cycles=1, energy_cap=1, power_cap=1, daa_price_vector=daa_price_vector)\n", - "step2_soc_ida, step2_cha_ida, step2_dis_ida, step2_cha_ida_close, step2_dis_ida_close, step2_profit_ida, step2_cha_daaida, step2_dis_daaida = optimizer_flex.step2_optimize_ida(n_cycles=1, energy_cap=1, power_cap=1, ida_price_vector=ida_price_vector, step1_cha_daa=step1_cha_daa, step1_dis_daa=step1_dis_daa)\n", - "step3_soc_idc, step3_cha_idc, step3_dis_idc, step3_cha_idc_close, step3_dis_idc_close, step3_profit_idc, step3_cha_daaidaidc, step3_dis_daaidaidc = optimizer_flex.step3_optimize_idc(n_cycles=1, energy_cap=1, power_cap=1, idc_price_vector=idc_price_vector, step2_cha_daaida=step2_cha_daaida, step2_dis_daaida=step2_dis_daaida)" + "step2_soc_ida, step2_cha_ida, step2_dis_ida, step2_cha_ida_close, step2_dis_ida_close, step2_profit_ida, step2_cha_daaida, step2_dis_daaida = optimizer_flex.step2_optimize_ida(n_cycles=1, energy_cap=1, power_cap=1, ida_price_vector=ida_price_vector)\n", + "step3_soc_idc, step3_cha_idc, step3_dis_idc, step3_cha_idc_close, step3_dis_idc_close, step3_profit_idc, step3_cha_daaidaidc, step3_dis_daaidaidc = optimizer_flex.step3_optimize_idc(n_cycles=1, energy_cap=1, power_cap=1, idc_price_vector=idc_price_vector)" ] }, { @@ -136,10 +136,11 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 12, "metadata": {}, "outputs": [], "source": [ + "# The revenue formula for each market are as fallows\n", "revenue_daa = np.sum(np.asarray(daa_price_vector) * (np.asarray(step1_dis_daa) - step1_cha_daa)) * 1/4\n", "revenue_ida = np.sum(np.asarray(ida_price_vector) * (np.asarray(step2_dis_ida) + step2_dis_ida_close - step2_cha_ida - step2_cha_ida_close)) * 1/4\n", "revenue_idc = np.sum(np.asarray(idc_price_vector) * (np.asarray(step3_dis_idc) + step3_dis_idc_close - step3_cha_idc - step3_cha_idc_close)) * 1/4\n", @@ -147,6 +148,19 @@ "revenue_total = revenue_daa+revenue_ida+revenue_idc" ] }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": {}, + "outputs": [], + "source": [ + "# The optimizer class contains revenue caculcation internally.\n", + "revenue_daa = optimizer_flex.revenue_daa()\n", + "revenue_ida = optimizer_flex.revenue_ida()\n", + "revenue_idc = optimizer_flex.revenue_idc()\n", + "revenue_total = optimizer_flex.revenue_total()" + ] + }, { "cell_type": "markdown", "metadata": {}, @@ -156,7 +170,7 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": 14, "metadata": {}, "outputs": [ { diff --git a/optimizer.py b/optimizer.py index efdf4c2..3ae2dee 100644 --- a/optimizer.py +++ b/optimizer.py @@ -177,6 +177,7 @@ def dis_daa_quarters_3_4_parity(model, q): # Define objective function and solve the optimization problem. # The objective is to maximize revenue from DA Auction trades over all possible charge-discharge schedules. + self.daa_price_vector = daa_price_vector model.obj = pyo.Objective(expr=sum(power_cap/4 * daa_price_vector[q-1] * ( model.dis_daa[q] - model.cha_daa[q]) for q in model.Q), sense=pyo.maximize) @@ -186,21 +187,21 @@ def dis_daa_quarters_3_4_parity(model, q): # Retrieve arrays of resulting optimal soc/charge/discharge schedules after the DA Auction: - step1_soc_daa = [model.soc[q].value for q in range( + self.step1_soc_daa = [model.soc[q].value for q in range( 1, len(daa_price_vector) + 1)] - step1_cha_daa = [model.cha_daa[q].value for q in range( + self.step1_cha_daa = [model.cha_daa[q].value for q in range( 1, len(daa_price_vector) + 1)] - step1_dis_daa = [model.dis_daa[q].value for q in range( + self.step1_dis_daa = [model.dis_daa[q].value for q in range( 1, len(daa_price_vector) + 1)] # Calculate profit from Day-ahead auction trades: - step1_profit_daa = sum([power_cap/4 * daa_price_vector[q] * ( - step1_dis_daa[q] - step1_cha_daa[q]) for q in range(len(daa_price_vector))]) + self.step1_profit_daa = sum([power_cap/4 * daa_price_vector[q] * ( + self.step1_dis_daa[q] - self.step1_cha_daa[q]) for q in range(len(daa_price_vector))]) - return(step1_soc_daa, step1_cha_daa, step1_dis_daa, step1_profit_daa) + return(self.step1_soc_daa, self.step1_cha_daa, self.step1_dis_daa, self.step1_profit_daa) - def step2_optimize_ida(self, n_cycles: int, energy_cap: int, power_cap: int, ida_price_vector: list, step1_cha_daa: list, step1_dis_daa: list): + def step2_optimize_ida(self, n_cycles: int, energy_cap: int, power_cap: int, ida_price_vector: list): """ Calculates optimal charge/discharge schedule on the intraday auction (ida) for a given 96-d ida_price_vector. @@ -292,43 +293,43 @@ def soc_step_constraint(model, q): """ The state of charge of each quarter equals the state if charge of the previous quarter plus charges minus discharges. (Constraint 2.5) """ - return model.soc[q+1] == model.soc[q] + power_cap/4 * (model.cha_ida[q] - model.dis_ida[q] + model.cha_ida_close[q] - model.dis_ida_close[q] + step1_cha_daa[q-1] - step1_dis_daa[q-1]) + return model.soc[q+1] == model.soc[q] + power_cap/4 * (model.cha_ida[q] - model.dis_ida[q] + model.cha_ida_close[q] - model.dis_ida_close[q] + self.step1_cha_daa[q-1] - self.step1_dis_daa[q-1]) def charge_cycle_limit(model): """ Sum of all charges has to be below the daily limit. (Constraint 2.6) """ - return ((np.sum(step1_cha_daa) + sum(model.cha_ida[q] for q in model.Q) - sum(model.dis_ida_close[q] for q in model.Q)) * power_cap/4 <= volume_limit) + return ((np.sum(self.step1_cha_daa) + sum(model.cha_ida[q] for q in model.Q) - sum(model.dis_ida_close[q] for q in model.Q)) * power_cap/4 <= volume_limit) def discharge_cycle_limit(model): """ Sum of all discharges has to be below the daily limit. (Constraint 2.7) """ - return ((np.sum(step1_dis_daa) + sum(model.dis_ida[q] for q in model.Q) - sum(model.cha_ida_close[q] for q in model.Q)) * power_cap/4 <= volume_limit) + return ((np.sum(self.step1_dis_daa) + sum(model.dis_ida[q] for q in model.Q) - sum(model.cha_ida_close[q] for q in model.Q)) * power_cap/4 <= volume_limit) def cha_close_logic(model, q): """ cha_ida_close can only close or reduce existing dis_daa positions. They can only be placed, where dis_daa positions exist. (Constraint 2.8) """ - return model.cha_ida_close[q] <= step1_dis_daa[q-1] + return model.cha_ida_close[q] <= self.step1_dis_daa[q-1] def dis_close_logic(model, q): """ dis_ida_close can only close or reduce existing cha_daa positions. They can only be placed, where cha_daa positions exist. (Constraint 2.9) """ - return model.dis_ida_close[q] <= step1_cha_daa[q-1] + return model.dis_ida_close[q] <= self.step1_cha_daa[q-1] def charge_rate_limit(model, q): """ Sum of cha_ida[q] and cha_daa[q] has to be less or equal to 1. (Constraint 2.10) """ - return model.cha_ida[q] + step1_cha_daa[q-1] <= 1 + return model.cha_ida[q] + self.step1_cha_daa[q-1] <= 1 def discharge_rate_limit(model, q): """ Sum of dis_ida[q] and dis_daa[q] has to be less or equal to 1. (Constraint 2.11) """ - return model.dis_ida[q] + step1_dis_daa[q-1] <= 1 + return model.dis_ida[q] + self.step1_dis_daa[q-1] <= 1 # Apply constraints on the model: @@ -353,6 +354,8 @@ def discharge_rate_limit(model, q): # Define objective function and solve the optimization problem # The objective is to maximize revenue from ID Auction trades over all possible charge-discharge schedules. + self.ida_price_vector = ida_price_vector + model.obj = pyo.Objective(expr=sum(ida_price_vector[q-1] * power_cap/4 * ( model.dis_ida[q] + model.dis_ida_close[q] - model.cha_ida[q] - model.cha_ida_close[q]) for q in model.Q), sense=pyo.maximize) @@ -361,32 +364,32 @@ def discharge_rate_limit(model, q): # Retrieve arrays of resulting optimal soc/charge/discharge schedules after the ID Auction: - step2_soc_ida = [model.soc[q].value for q in range( + self.step2_soc_ida = [model.soc[q].value for q in range( 1, len(ida_price_vector) + 1)] - step2_cha_ida = [model.cha_ida[q].value for q in range( + self.step2_cha_ida = [model.cha_ida[q].value for q in range( 1, len(ida_price_vector) + 1)] - step2_dis_ida = [model.dis_ida[q].value for q in range( + self.step2_dis_ida = [model.dis_ida[q].value for q in range( 1, len(ida_price_vector) + 1)] - step2_cha_ida_close = [model.cha_ida_close[q].value for q in range( + self.step2_cha_ida_close = [model.cha_ida_close[q].value for q in range( 1, len(ida_price_vector) + 1)] - step2_dis_ida_close = [model.dis_ida_close[q].value for q in range( + self.step2_dis_ida_close = [model.dis_ida_close[q].value for q in range( 1, len(ida_price_vector) + 1)] # Calculate profit from Day-ahead auction trades: - step2_profit_ida = np.sum(((np.asarray(step2_dis_ida) + step2_dis_ida_close) - ( - np.asarray(step2_cha_ida) + step2_cha_ida_close)) * ida_price_vector) * power_cap/4 + self.step2_profit_ida = np.sum(((np.asarray(self.step2_dis_ida) + self.step2_dis_ida_close) - ( + np.asarray(self.step2_cha_ida) + self.step2_cha_ida_close)) * ida_price_vector) * power_cap/4 # Calculate total physical charge discharge schedules of combined day-ahead and intraday auction trades: - step2_cha_daaida = np.asarray( - step1_cha_daa) - step2_dis_ida_close + step2_cha_ida - step2_dis_daaida = np.asarray( - step1_dis_daa) - step2_cha_ida_close + step2_dis_ida + self.step2_cha_daaida = np.asarray( + self.step1_cha_daa) - self.step2_dis_ida_close + self.step2_cha_ida + self.step2_dis_daaida = np.asarray( + self.step1_dis_daa) - self.step2_cha_ida_close + self.step2_dis_ida - return(step2_soc_ida, step2_cha_ida, step2_dis_ida, step2_cha_ida_close, step2_dis_ida_close, step2_profit_ida, step2_cha_daaida, step2_dis_daaida) + return(self.step2_soc_ida, self.step2_cha_ida, self.step2_dis_ida, self.step2_cha_ida_close, self.step2_dis_ida_close, self.step2_profit_ida, self.step2_cha_daaida, self.step2_dis_daaida) - def step3_optimize_idc(self, n_cycles: int, energy_cap: int, power_cap: int, idc_price_vector: list, step2_cha_daaida: list, step2_dis_daaida: list): + def step3_optimize_idc(self, n_cycles: int, energy_cap: int, power_cap: int, idc_price_vector: list): """ Calculates optimal charge/discharge schedule on the intraday continuous (idc) for a given 96-d idc_price_vector. @@ -478,43 +481,43 @@ def soc_step_constraint(model, q): """ The state of charge of each quarter equals the state if charge of the previous quarter plus charges minus discharges. (Constraint 3.5) """ - return model.soc[q+1] == model.soc[q] + power_cap/4 * (model.cha_idc[q] - model.dis_idc[q] + model.cha_idc_close[q] - model.dis_idc_close[q] + step2_cha_daaida[q-1] - step2_dis_daaida[q-1]) + return model.soc[q+1] == model.soc[q] + power_cap/4 * (model.cha_idc[q] - model.dis_idc[q] + model.cha_idc_close[q] - model.dis_idc_close[q] + self.step2_cha_daaida[q-1] - self.step2_dis_daaida[q-1]) def charge_cycle_limit(model): """ Sum of all charges has to be below the daily limit. (Constraint 3.6) """ - return (np.sum(step2_dis_daaida) + sum(model.dis_idc[q] for q in model.Q) - sum(model.cha_idc_close[q] for q in model.Q)) * power_cap/4 <= volume_limit + return (np.sum(self.step2_dis_daaida) + sum(model.dis_idc[q] for q in model.Q) - sum(model.cha_idc_close[q] for q in model.Q)) * power_cap/4 <= volume_limit def discharge_cycle_limit(model): """ Sum of all discharges has to be below the daily limit. (Constraint 3.7) """ - return (np.sum(step2_cha_daaida) + sum(model.cha_idc[q] for q in model.Q) - sum(model.dis_idc_close[q] for q in model.Q)) * power_cap/4 <= volume_limit + return (np.sum(self.step2_cha_daaida) + sum(model.cha_idc[q] for q in model.Q) - sum(model.dis_idc_close[q] for q in model.Q)) * power_cap/4 <= volume_limit def cha_close_logic(model, q): """ cha_idc_close can only close or reduce existing dis_daaida positions. They can only be placed, where dis_daaida positions exist. (Constraint 3.8) """ - return model.cha_idc_close[q] <= step2_dis_daaida[q-1] + return model.cha_idc_close[q] <= self.step2_dis_daaida[q-1] def dis_close_logic(model, q): """ dis_idc_close can only close or reduce existing cha_daaida positions. They can only be placed, where cha_daaida positions exist. (Constraint 3.9) """ - return model.dis_idc_close[q] <= step2_cha_daaida[q-1] + return model.dis_idc_close[q] <= self.step2_cha_daaida[q-1] def charge_rate_limit(model, q): """ Sum of cha_idc[q] and cha_daaida[q] has to be less or equal to 1. (Constraint 3.10) """ - return model.cha_idc[q] + step2_cha_daaida[q-1] <= 1 + return model.cha_idc[q] + self.step2_cha_daaida[q-1] <= 1 def discharge_rate_limit(model, q): """ Sum of dis_idc[q] and dis_daaida[q] has to be less or equal to 1. (Constraint 3.11) """ - return model.dis_idc[q] + step2_dis_daaida[q-1] <= 1 + return model.dis_idc[q] + self.step2_dis_daaida[q-1] <= 1 # Apply constraints on the model: @@ -540,6 +543,8 @@ def discharge_rate_limit(model, q): # Define objective function and solve the optimization problem # The objective is to maximize revenue from ID Continuous trades over all possible charge-discharge schedules. + self.idc_price_vector = idc_price_vector + model.obj = pyo.Objective(expr=sum([idc_price_vector[q-1] * power_cap/4 * ( model.dis_idc[q]+model.dis_idc_close[q]-model.cha_idc[q]-model.cha_idc_close[q]) for q in model.Q]), sense=pyo.maximize) @@ -548,27 +553,50 @@ def discharge_rate_limit(model, q): # Retrieve arrays of resulting optimal soc/charge/discharge schedules after the ID Auction: - step3_soc_idc = [model.soc[q].value for q in range( + self.step3_soc_idc = [model.soc[q].value for q in range( 1, len(idc_price_vector) + 1)] - step3_cha_idc = [model.cha_idc[q].value for q in range( + self.step3_cha_idc = [model.cha_idc[q].value for q in range( 1, len(idc_price_vector) + 1)] - step3_dis_idc = [model.dis_idc[q].value for q in range( + self.step3_dis_idc = [model.dis_idc[q].value for q in range( 1, len(idc_price_vector) + 1)] - step3_cha_idc_close = [model.cha_idc_close[q].value for q in range( + self.step3_cha_idc_close = [model.cha_idc_close[q].value for q in range( 1, len(idc_price_vector) + 1)] - step3_dis_idc_close = [model.dis_idc_close[q].value for q in range( + self.step3_dis_idc_close = [model.dis_idc_close[q].value for q in range( 1, len(idc_price_vector) + 1)] # Calculate profit from Day-ahead auction trades: - step3_profit_idc = np.sum(((np.asarray(step3_dis_idc) + step3_dis_idc_close) - ( - np.asarray(step3_cha_idc) + step3_cha_idc_close)) * idc_price_vector) * power_cap/4 + self.step3_profit_idc = np.sum(((np.asarray(self.step3_dis_idc) + self.step3_dis_idc_close) - ( + np.asarray(self.step3_cha_idc) + self.step3_cha_idc_close)) * idc_price_vector) * power_cap/4 # Calculate total physical charge discharge schedules of combined day-ahead and intraday auction trades: - step3_cha_daaidaidc = np.asarray( - step2_cha_daaida) - step3_dis_idc_close + step3_cha_idc - step3_dis_daaidaidc = np.asarray( - step2_dis_daaida) - step3_cha_idc_close + step3_dis_idc + self.step3_cha_daaidaidc = np.asarray( + self.step2_cha_daaida) - self.step3_dis_idc_close + self.step3_cha_idc + self.step3_dis_daaidaidc = np.asarray( + self.step2_dis_daaida) - self.step3_cha_idc_close + self.step3_dis_idc + + return(self.step3_soc_idc, self.step3_cha_idc, self.step3_dis_idc, self.step3_cha_idc_close, self.step3_dis_idc_close, self.step3_profit_idc, self.step3_cha_daaidaidc, self.step3_dis_daaidaidc) + + # Calculates the revenue of DAA, IDA, and IDC market. And total revenue. + # Each revenue formala can be found in https://github.com/FlexPwr/bess-optimizer/blob/main/mathematical_formulation.pdf + + def revenue_daa(self) -> float: + # Calculate the revenue of DAA market + + return(np.sum(np.asarray(self.daa_price_vector) * (np.asarray(self.step1_dis_daa) - self.step1_cha_daa)) * 1/4) + + def revenue_ida(self) -> float: + # Calculate the revenue of IDA market + + return(np.sum(np.asarray(self.ida_price_vector) * (np.asarray(self.step2_dis_ida) + self.step2_dis_ida_close - self.step2_cha_ida - self.step2_cha_ida_close)) * 1/4) + + def revenue_idc(self) -> float: + # Calculate the revenue of IDC market + + return(np.sum(np.asarray(self.idc_price_vector) * (np.asarray(self.step3_dis_idc) + self.step3_dis_idc_close - self.step3_cha_idc - self.step3_cha_idc_close)) * 1/4) + + def revenue_total(self) -> float: + # Calculate the total revenue - return(step3_soc_idc, step3_cha_idc, step3_dis_idc, step3_cha_idc_close, step3_dis_idc_close, step3_profit_idc, step3_cha_daaidaidc, step3_dis_daaidaidc) + return(self.revenue_daa() + self.revenue_ida() + self.revenue_idc())