fix: over-abstraction flow for RunModel.Supervisor

Refs #144

R/RunModel.InputsModel.R

@@ -60,12 +60,14 @@ RunModel.InputsModel <- function(x = NULL,
     # Upstream basins and Reservoir are launch directly
     OutputsModel <- FUN_MOD(x, RunOptions, Param)
   } else {
-    # Intermediate basins (other than reservoir) are launch with SD capabilities
+    # Intermediate basins (other than reservoir) are launched with SD capabilities
     if (!is.null(x$UpstreamNodes) & !inherits(x, "SD")) {
       # For calibration of node with diversion
       class(x) <- c(class(x), "SD")
     }
     OutputsModel <- airGR::RunModel(x, RunOptions, Param, FUN_MOD)
+    OutputsModel <- calcOverAbstraction(OutputsModel, FALSE)
+    OutputsModel$RunOptions <- calcOverAbstraction(OutputsModel$RunOptions, TRUE)
   }
   OutputsModel$RunOptions$TimeStep <- RunOptions$FeatFUN_MOD$TimeStep
   if (is.null(OutputsModel$Qsim_m3)) {
@@ -80,8 +82,6 @@ RunModel.InputsModel <- function(x = NULL,
   if (x$hasDiversion && !x$isReservoir) {
     OutputsModel <- RunModel_Diversion(x, RunOptions, OutputsModel)
   }
-  OutputsModel <- calcOverAbstraction(OutputsModel, FALSE)
-  OutputsModel$RunOptions <- calcOverAbstraction(OutputsModel$RunOptions, TRUE)
   return(OutputsModel)
 }
 
@@ -144,29 +144,3 @@ calc_Qdiv<- function(Qnat, Qdiv, Qmin) {
   }
   return(list(Qsim = Qsim, Qdiv = Qnat - Qsim))
 }
-
-
-#' Cap negative `OutputsModel$Qsim_m3` to zero and fill `OutputsModel$Qover_m3`
-#' with over-abstracted volumes
-#'
-#' @param O Either `OutputsModel` or `OutputsModel$RunOptions` (for warm-up Qsim)
-#' @param WarmUp `TRUE` if `O` is `OutputsModel$RunOptions`
-#'
-#' @return Modified `OutputsModel` or `OutputsModel$RunOptions`
-#' @noRd
-#'
-calcOverAbstraction <- function(O, WarmUp) {
-  f <- list(sim = "Qsim_m3", over = "Qover_m3")
-  if(WarmUp) {
-    f <- lapply(f, function(x) paste0("WarmUp", x))
-  }
-  if (!is.null(O[[f$sim]])) {
-    if (any(!is.na(O[[f$sim]]) & O[[f$sim]] < 0)) {
-      O[[f$over]] <- rep(0, length(O[[f$sim]]))
-      O[[f$over]][O[[f$sim]] < 0] <- - O[[f$sim]][!is.na(O[[f$sim]]) & O[[f$sim]] < 0]
-      O[[f$sim]][!is.na(O[[f$sim]]) & O[[f$sim]] < 0] <- 0
-    }
-  }
-  return(O)
-}
-

R/RunModel.SD.R

-#' Run a semi-distributed model from rainfall-runoff model outputs
+#' Run a hydraulic routing model from rainfall-runoff model outputs
 #'
 #' @inheritParams airGR::RunModel_Lag
 #' @param x \[object of class `InputsModel`\] used as `InputsModel` parameter for [airGR::RunModel_Lag]
@@ -17,6 +17,8 @@ RunModel.SD <- function(x, RunOptions, Param, QcontribDown, ...) {
                                         RunOptions = RunOptions,
                                         Param = Param[1],
                                         QcontribDown = QcontribDown)
+    OutputsModel <- calcOverAbstraction(OutputsModel, FALSE)
+    OutputsModel$RunOptions <- calcOverAbstraction(OutputsModel$RunOptions, TRUE)
   }
   OutputsModel$RunOptions$TimeStep <- RunOptions$FeatFUN_MOD$TimeStep
   return(OutputsModel)

R/RunModel.Supervisor.R

@@ -48,9 +48,6 @@ RunModel.Supervisor <- function(x, RunOptions, Param, ...) {
     updateQupstream.Supervisor(x, id, IndPeriod_Run)
   }
 
-  # Store OutputsModel for step by step simulation
-  x$storedOutputs <- initStoredOutputs(x)
-
   # Initialization of model states by running the model with no supervision on warm-up period
   RunOptionsWarmUp <- RunOptions
   for(id in names(x$InputsModel)) {
@@ -75,14 +72,17 @@ RunModel.Supervisor <- function(x, RunOptions, Param, ...) {
 
   # Set Outputs to archive for final restitution
   outputVars <- lapply(SD_Ids, function(id) {
-    ov <- "Qsim_m3"
+    ov <- c("Qsim_m3", "Qover_m3")
     if (x$InputsModel[[id]]$hasDiversion) {
       ov <- c(ov, "Qdiv_m3", "Qnat")
-    } else if (x$InputsModel[[id]]$isReservoir) {
+    }
+    if (x$InputsModel[[id]]$isReservoir) {
       ov <- c(ov, "Qinflows_m3", "Vsim")
     }
     return(ov)
   })
+  # Store OutputsModel for step by step simulation
+  x$storedOutputs <- initStoredOutputs(x, outputVars)
 
   message("Processing: 0%", appendLF = FALSE)
   iProgressSteps <- round(length(lSuperTS) * seq(0.1, 0.9, 0.1))
@@ -152,6 +152,7 @@ RunModel.Supervisor <- function(x, RunOptions, Param, ...) {
   return(x$OutputsModel)
 }
 
+
 updateQupstream.Supervisor <- function(x, id, iTS) {
   downId <- x$InputsModel[[id]]$down
   if(!is.null(x$InputsModel[[downId]])) {

R/utils.RunModel.R

@@ -64,3 +64,27 @@ serializeIniStates <- function(IniStates) {
   unlist(IniStates)
 }
 
+
+#' Cap negative `OutputsModel$Qsim_m3` to zero and fill `OutputsModel$Qover_m3`
+#' with over-abstracted volumes
+#'
+#' @param O Either `OutputsModel` or `OutputsModel$RunOptions` (for warm-up Qsim)
+#' @param WarmUp `TRUE` if `O` is `OutputsModel$RunOptions`
+#'
+#' @return Modified `OutputsModel` or `OutputsModel$RunOptions`
+#' @noRd
+#'
+calcOverAbstraction <- function(O, WarmUp) {
+  f <- list(sim = "Qsim_m3", over = "Qover_m3")
+  if(WarmUp) {
+    f <- lapply(f, function(x) paste0("WarmUp", x))
+  }
+  if (!is.null(O[[f$sim]])) {
+    O[[f$over]] <- rep(0, length(O[[f$sim]]))
+    if (any(!is.na(O[[f$sim]]) & O[[f$sim]] < 0)) {
+      O[[f$over]][O[[f$sim]] < 0] <- - O[[f$sim]][!is.na(O[[f$sim]]) & O[[f$sim]] < 0]
+      O[[f$sim]][!is.na(O[[f$sim]]) & O[[f$sim]] < 0] <- 0
+    }
+  }
+  return(O)
+}

R/utils.Supervisor.R

@@ -97,30 +97,21 @@ doSupervision <- function(supervisor) {
 }
 
 
-initStoredOutputs <- function(x) {
-  np <- getAllNodesProperties(x$griwrm)
-  so <- list()
-  so$QcontribDown <- do.call(
+initStoredOutputs <- function(x, outputVars) {
+  QcontribDown <- do.call(
     cbind,
     lapply(x$OutputsModel, "[[", "Qsim")
   )
-  so$Qsim_m3 <- do.call(
-    cbind,
-    lapply(x$OutputsModel, "[[", "Qsim_m3")
-  )
-  if (sum(np$Diversion) > 0) {
-    # Outputs of Diversion nodes
-    so$Qdiv_m3 <- so$Qsim_m3[, np$id[np$Diversion], drop = FALSE] * NA
-    so$Qnat <- so$Qdiv_m3
-  }
-  if (sum(np$Reservoir) > 0) {
-    # Specific Outputs of RunModel_Reservoir
-    so$Vsim <- matrix(rep(NA, sum(np$Reservoir) * nrow(so$Qsim_m3)),
-                      nrow = nrow(so$Qsim_m3))
-    colnames(so$Vsim) <- np$id[np$Reservoir]
-    so$Qinflows_m3 <- so$Vsim
-    # Add columns Qsim_m3 at reservoir (out of the scope of GR models calculated above)
-    so$Qsim_m3 <- cbind(so$Qsim_m3, so$Vsim)
-  }
+  so <- lapply(setNames(nm = unique(unlist(outputVars))), function(ov) {
+    s <- sapply(outputVars, function(y) "Qsim_m3" %in% y)
+    ids <- names(s)[s]
+    if (length(ids) > 0) {
+      m <- matrix(NA, nrow = nrow(QcontribDown), ncol = length(ids))
+      colnames(m) <- ids
+      return(m)
+    }
+    return(NULL)
+  })
+  so$QcontribDown <- QcontribDown
   return(so)
 }

tests/testthat/test-RunModel.R

@@ -261,4 +261,8 @@ test_that("RunModel should return water deficit (Qover_m3)", {
   for(x in ls(e)) assign(x, get(x, e))
   expect_false(any(OM_GriwrmInputs$`54001`$Qsim_m3 < 0))
   expect_true(all(OM_GriwrmInputs$`54001`$Qover_m3 >= 0))
+  sv <- CreateSupervisor(InputsModel)
+  OM_sv <- RunModel(sv, RunOptions, ParamMichel)
+  expect_equal(OM_sv$`54001`$Qsim_m3, OM_GriwrmInputs$`54001`$Qsim_m3)
+  expect_equal(OM_sv$`54001`$Qover_m3, OM_GriwrmInputs$`54001`$Qover_m3)
 })