New function getBellmanValuesWithPlaia (#20)

* New function getBellmanValuesWithPlaia

* Avoid errors

* Cvar in getBellmanValuesWithPlaia

* Update function get bellman plaia

* Fix errors

* Refactor plaia (#23)

* validate_and_normalize_areas

* validate_api_study_for_plaia

* prepare_areas_for_simulation

* run_plaia_simulation and download_output_zip

* write_api_file

* extract_from_zip

* Update dowload_output_zip

* Fix check

* fix

* prepare_and_launch_plaia

* Fux error

* Fix erros zipfile

* other_options in prepare_and_launch_plaia

* Fix bellman values

* Improve bellman overflow (#24)

* Simplify bellman values interpolation

* Simplify next states values

* Reorganize build_all_possible_decisions

* Better overflow

* Improve documentation (#25)

* Rewrite README

* New table for methods overview

* Remove section multistock

* Overflow comment

* Increment version number to 2.2.0

* Argument final level for all functions

* New function to wait for task end

* Reduce time

* Update NEWS.md

Author: Juliette-Gerbaux

DESCRIPTION

@@ -1,7 +1,7 @@
 Type: Package
 Package: antaresWaterValues
 Title: Generating water values for Antares
-Version: 2.1.1
+Version: 2.2.0
 Authors@R: c(
     person("Dhia", "Gharsallaoui", role = "aut"),
     person("Juliette", "Gerbaux", email = "juliette.gerbaux@rte-france.com",role = c("aut","cre")),
@@ -42,6 +42,7 @@ Suggests:
     spsComps (>= 0.3.0),
     spsUtil,
     testthat (>= 3.0.0),
+    yaml,
     zoo
 VignetteBuilder: 
     knitr

NAMESPACE

@@ -7,6 +7,7 @@ export(changeHydroManagement)
 export(constraint_generator)
 export(getBellmanValuesFromOneSimulationMultistock)
 export(getBellmanValuesSequentialMultiStockWithPlaia)
+export(getBellmanValuesWithPlaia)
 export(getPumpEfficiency)
 export(get_Reward)
 export(get_inflow)

NEWS.md

@@ -1,3 +1,15 @@
+# antaresWaterValues 2.2.0
+
+## New features
+
+* New function to compute Bellman values directly with plaia. 
+* Verification of task ending in API mode.
+* Argument final_level on almost all Bellman values functions.
+
+## Breaking change
+
+* New way to handle overflow in dynamic programming. It can be done anytime but is limited by weekly inflows and is penalized by spillage cost. 
+
 # antaresWaterValues 2.1.0
 
 ## New features

R/Bellman.R

@@ -47,12 +47,9 @@
     # Get interpolation function of rewards for each possible transition for each MC year
     f_reward_year <- get_reward_interpolation(Data_week)
 
-    #Get interpolation function of next Bellman values
-    f_next_value <- get_bellman_values_interpolation(next_state,next_week_values_l,mcyears)
-
     # Build a data.table from Data_week that list for each state and each MC year, the possible transitions
-    df_SDP <- build_all_possible_decisions(Data_week,decision_space,f_next_value,
-                                           mcyears,lvl_high,lvl_low,E_max,P_max,
+    df_SDP <- build_all_possible_decisions(Data_week,decision_space,
+                                           lvl_high,lvl_low,E_max,P_max,
                                            next_week_values_l,niveau_max,overflow_cost,
                                            next_states = next_state)
 

R/Bellman_subfunctions.R

@@ -15,30 +15,6 @@ get_reward_interpolation <- function(Data_week){
   return(f_reward_year)
 }
 
-#' Create approximation of Bellman function for next week for each scenario, used in \code{Bellman} and
-#' in \code{getOptimalTrend}
-#'
-#' @param next_state possible states for next week
-#' @param next_week_values Vector of Bellman values to interpolate
-#' @param mcyears Vector of monte carlo years to take into account
-#'
-#' @return List of \code{stats::approxfun} for each scenario
-#' @keywords internal
-get_bellman_values_interpolation <- function(next_state,next_week_values,mcyears){
-
-  df_next_week <- data.frame(years = mcyears,
-                             next_state = next_state,
-                             next_value = next_week_values)
-
-  f_next_value <- c()
-  for (year in mcyears){
-    df <- dplyr::filter(df_next_week, .data$years==year)
-    f <- stats::approxfun(df$next_state, df$next_value)#,yleft = 0, yright=0)
-    f_next_value <- c(f_next_value,f)
-  }
-
-  return(f_next_value)
-}
 
 #' Create a data frame with all possible transition and associated next state for
 #' each scenario and each state for one particular week, used in \code{Bellman}.
@@ -49,9 +25,6 @@ get_bellman_values_interpolation <- function(next_state,next_week_values,mcyears
 #' as weel as inflows
 #' @param decision_space Simulation constraints values for the week
 #'  generated by \code{get_Reward} or \code{runWaterValuesSimulation}
-#' @param f_next_value List of approximated Bellman functions for the next week,
-#' generated by \code{get_bellman_values_interpolation}
-#' @param mcyears Vector of monte carlo years to take into account
 #' @param lvl_high Double, top rule curve for the end of the week in MWh
 #' @param lvl_low Double, bottom rule curve for the end of the week in MWh
 #' @param E_max Double, weekly maximum turbining power for the week
@@ -65,42 +38,47 @@ get_bellman_values_interpolation <- function(next_state,next_week_values,mcyears
 #' next_state, control (transition to next state) and
 #' next value (Bellman value of next state)
 #' @keywords internal
-build_all_possible_decisions <- function(Data_week,decision_space,f_next_value,
-                                         mcyears,lvl_high,lvl_low,E_max,P_max,
+build_all_possible_decisions <- function(Data_week,decision_space,
+                                         lvl_high,lvl_low,E_max,P_max,
                                          next_week_values,niveau_max,overflow_cost,
                                          next_states){
 
-  df_next_week <- data.frame(years = Data_week$years,
-                             next_state = next_states,
-                             next_value = next_week_values)
+  df_next_week <- dplyr::distinct(data.frame(
+    next_state = next_states,
+    next_value = next_week_values
+  ))
 
-  future_states <- Data_week %>%
-    dplyr::inner_join(df_next_week,by="years", relationship="many-to-many") %>%
-    dplyr::mutate(control = -.data$next_state+.data$states+.data$hydroStorage)
+  #Get interpolation function of next Bellman values
+  f_next_value <- stats::approxfun(df_next_week$next_state, df_next_week$next_value)
 
-  control_possible <- Data_week  %>%
-    dplyr::right_join(decision_space,by=c("years"="mcYear"), relationship="many-to-many") %>%
-    dplyr::rename("control"="u") %>%
-    dplyr::mutate(next_state=dplyr::if_else(.data$states+.data$hydroStorage-.data$control>niveau_max,niveau_max,
-                                            .data$states+.data$hydroStorage-.data$control)) %>%
-    dplyr::mutate(overflow=dplyr::if_else(.data$states+.data$hydroStorage-.data$control>niveau_max,
-                                          .data$states+.data$hydroStorage-.data$control-niveau_max,0)) %>%
-    dplyr::mutate(next_value=mapply(function(y,x)f_next_value[[which(y==mcyears)]](x), .data$years, .data$next_state)+
-                    .data$overflow*overflow_cost) %>%
-    dplyr::select(-c("overflow"))
-
-  control_min <- Data_week %>%
-    dplyr::mutate(next_state=lvl_high) %>%
-    dplyr::mutate(control = -.data$next_state+.data$states+.data$hydroStorage) %>%
-    dplyr::mutate(next_value=mapply(function(y,x)f_next_value[[which(y==mcyears)]](x), .data$years, .data$next_state))
+  possible_states = data.frame(next_state = c(df_next_week$next_state,
+                                              lvl_low,lvl_high))
 
-  control_max <- Data_week %>%
-    dplyr::mutate(next_state=lvl_low) %>%
+  possible_control <- Data_week %>%
+    dplyr::cross_join(possible_states) %>%
     dplyr::mutate(control = -.data$next_state+.data$states+.data$hydroStorage) %>%
-    dplyr::mutate(next_value=mapply(function(y,x)f_next_value[[which(y==mcyears)]](x), .data$years, .data$next_state))
+    dplyr::filter((-P_max<=.data$control)&(.data$control<=E_max))
 
-  df_SDP <- dplyr::bind_rows(future_states, control_possible, control_min, control_max) %>%
-    dplyr::filter((-P_max<=.data$control)&(.data$control<=E_max)&(.data$next_state>=0))
+  possible_control <- Data_week  %>%
+    dplyr::right_join(decision_space,by=c("years"="mcYear"), relationship="many-to-many") %>%
+    dplyr::rename("control"="u") %>%
+    dplyr::mutate(next_state = -.data$control+.data$states+.data$hydroStorage) %>%
+    rbind(possible_control)
+
+  df_SDP = possible_control %>%
+    dplyr::select(-c("next_state")) %>%
+    dplyr::cross_join(possible_states) %>%
+    dplyr::mutate(overflow = -.data$next_state+.data$states+.data$hydroStorage - .data$control) %>%
+    dplyr::filter(0<=.data$overflow,.data$overflow<=.data$hydroStorage) %>%
+    rbind(dplyr::mutate(possible_control, overflow = 0)) %>%
+    rbind(dplyr::mutate(possible_control,
+                        overflow = .data$hydroStorage,
+                        next_state = .data$next_state-.data$overflow))
+
+  df_SDP <- df_SDP %>%
+    dplyr::filter(.data$next_state>=0, .data$next_state<=niveau_max)%>%
+    dplyr::mutate(next_value = f_next_value(.data$next_state) + .data$overflow*overflow_cost) %>%
+    dplyr::select(-c("overflow"))
 
   return(df_SDP)
 }

R/Grid_Matrix.R

@@ -1,7 +1,8 @@
 #' Compute Bellman values
 #'
 #' Compute Bellman values and water values with dynamic programming based on reward functions computed with \code{get_Reward()}.
-#' Mode information about this function in \code{vignette("grid_Matrix-parameters")}.
+#' Mode information about this function in \code{vignette("grid_Matrix-parameters")}. In the backward recursion, overflow is accepted but it
+#' can't exceed weekly inflows and is penalized with the spillage cost.
 #'
 #' @inheritParams runWaterValuesSimulation
 #' @param reward_db Reward functions generated by the function \code{get_Reward()}.
@@ -15,7 +16,7 @@
 #' @param max_hydro_weekly Maximum weekly pumping and generating power generated by the function \code{get_max_hydro()} with \code{timeStep="weekly"}.
 #' @param force_final_level Binary. Whether final level should be constrained.
 #' @param final_level Double. Final level (in percent between 0 and 100) if final level is constrained. If you want initial level, use \code{get_initial_level()}.
-#' @param penalty_final_level_low Double. Penalties for both bottom rule curve to constrain final level.
+#' @param penalty_final_level_low Double. Penalties for bottom rule curve to constrain final level.
 #' @param penalty_final_level_high Double. Penalties for top rule curve to constrain final level.
 #' @param until_convergence Binary. \code{TRUE} to repeat cycles until convergence of water values or
 #'  attending the limit.

R/iterations_simulation_DP.R

@@ -18,6 +18,7 @@
 #' between the reservoir capacity and zero
 #' @param cvar_value from 0 to 1. the probability used in cvar method
 #' @param penalty_final_level Penalties (for both bottom and top rule curves) to constrain final level
+#' @param final_level Double. Final level (in percent between 0 and 100). Initial level computed by \code{get_initial_level()} by default.
 #' @param df_previous_cut Data frame containing previous estimations of cuts
 #'
 #' @export
@@ -28,6 +29,7 @@ calculateBellmanWithIterativeSimulations <- function(area,pumping, pump_eff=1,op
                                                      path_solver,
                                                      states_step_ratio=1/50,
                                                      cvar_value = 1,
+                                                     final_level = NULL,
                                                      penalty_final_level = NULL,
                                                      df_previous_cut = NULL){
 
@@ -76,7 +78,9 @@ calculateBellmanWithIterativeSimulations <- function(area,pumping, pump_eff=1,op
   controls <- tidyr::drop_na(controls) %>%
     dplyr::cross_join(data.frame(mcYear=mcyears))
 
-  final_level = get_initial_level(area,opts)
+  if (is.null(final_level)){
+    final_level = get_initial_level(area,opts)
+  }
 
   if (!is.null(df_previous_cut)){
     controls = df_previous_cut %>%
@@ -512,6 +516,7 @@ getNewConstraint <- function(level_init,watervalues,mcyears,
 #' less to it. (quantile(cvar_value))
 #' @param penalty_final_level Penalties (for both bottom and top rule curves) to constrain final level
 #' @param initial_traj Initial trajectory (used for other storages)
+#' @param list_final_level List of double. For each storage, final level (in percent between 0 and 100) if final level is constrained. Initial level computed by \code{get_initial_level()} by default.
 #' @param list_areas_to_compute Vector of character. Areas for which to compute Bellman values. If \code{NULL}, all areas in \code{list_areas} are used.
 #' @inheritParams runWaterValuesSimulationMultiStock
 #' @inheritParams calculateBellmanWithIterativeSimulations
@@ -524,6 +529,7 @@ calculateBellmanWithIterativeSimulationsMultiStock <- function(list_areas,list_p
                                                                path_solver,
                                                                states_step_ratio=1/50,
                                                                cvar_value = 1,
+                                                               list_final_level = NULL,
                                                                penalty_final_level = NULL,
                                                                initial_traj = NULL,
                                                                df_previous_cut = NULL,
@@ -549,14 +555,13 @@ calculateBellmanWithIterativeSimulationsMultiStock <- function(list_areas,list_p
   assertthat::assert_that(opts$antaresVersion>=870,
                             msg = "Scenarization of rhs of binding constraints not available with the version of Antares. Update the study to 8.7.0 or don't scenarize control values.")
 
-  list_areas = tolower(list_areas)
+  validate_and_normalize_areas(list_areas,opts)
 
   list_backup = list()
   list_inflow = list()
   list_capacity = list()
   for (j in seq_along(list_areas)){
     area = list_areas[[j]]
-    check_area_name(area = area, opts = opts)
 
     list_backup[[j]] = getBackupData(area,mcyears,opts)
 
@@ -620,7 +625,12 @@ calculateBellmanWithIterativeSimulationsMultiStock <- function(list_areas,list_p
     a = area
     pumping <- list_pumping[area]
     pump_eff <- list_efficiency[area]
-    final_level <- get_initial_level(area,opts)
+    if (is.null(list_final_level)){
+      final_level <- get_initial_level(area,opts)
+    } else {
+      final_level <- list_final_level[[area]]
+    }
+
 
     max_hydro <- get_max_hydro(area, opts)
     max_hydro_weekly <- max_hydro %>%

R/multistock_1simulation.R

@@ -7,6 +7,7 @@
 #'
 #' @inheritParams runWaterValuesSimulationMultiStock
 #' @inheritParams Grid_Matrix
+#' @param list_final_level List of double. For each storage, final level (in percent between 0 and 100) if final level is constrained. Initial level computed by \code{get_initial_level()} by default.
 #' @param write_vu Binary. True to write water values in the Antares study.
 #' @param simu List of simulation parameters returned by the function \code{antaresRead::setSimulationPath()} with the simulation selected from which to use the storage trajectory to run the simulation.
 #'
@@ -19,6 +20,7 @@ getBellmanValuesFromOneSimulationMultistock <- function(opts,
                                                         list_areas,
                                                         list_pumping,
                                                         list_efficiency,
+                                                        list_final_level = NULL,
                                                         force_final_level,
                                                         penalty_final_level_low,
                                                         penalty_final_level_high,
@@ -92,7 +94,11 @@ getBellmanValuesFromOneSimulationMultistock <- function(opts,
   for (a in list_areas) {
     pumping <- list_pumping[a]
     pump_eff <- list_efficiency[a]
-    final_level <- get_initial_level(a, opts)
+    if (!is.null(list_final_level)){
+      final_level <- list_final_level[[a]]
+    } else {
+      final_level <- get_initial_level(a, opts)
+    }
 
     print(a)
     local_sim_values <- simulation_res$simulation_values %>%