From dccbb2da79466f5db7dabb380c2889a0f3274eab Mon Sep 17 00:00:00 2001
From: David <david.dorchies@inrae.fr>
Date: Fri, 29 Mar 2024 18:05:41 +0100
Subject: [PATCH] feat: Computation of Vsim and Qdiv for Diversion in Reservoir
Refs #146
---
R/RunModel.InputsModel.R | 2 +-
R/RunModel_Reservoir.R | 13 ++++
tests/testthat/helper_1_RunModel.R | 87 ++++++++++++------------
tests/testthat/test-RunModel_Reservoir.R | 15 +++-
4 files changed, 69 insertions(+), 48 deletions(-)
diff --git a/R/RunModel.InputsModel.R b/R/RunModel.InputsModel.R
index f36e4b0..8fb7a70 100644
--- a/R/RunModel.InputsModel.R
+++ b/R/RunModel.InputsModel.R
@@ -70,7 +70,7 @@ RunModel.InputsModel <- function(x = NULL,
OutputsModel$RunOptions$WarmUpQsim_m3 <-
OutputsModel$RunOptions$WarmUpQsim * sum(x$BasinAreas, na.rm = TRUE) * 1e3
}
- if (x$hasDiversion) {
+ if (x$hasDiversion && !x$isReservoir) {
OutputsModel <- RunModel_Diversion(x, RunOptions, OutputsModel)
}
return(OutputsModel)
diff --git a/R/RunModel_Reservoir.R b/R/RunModel_Reservoir.R
index 0fff58a..3e20530 100644
--- a/R/RunModel_Reservoir.R
+++ b/R/RunModel_Reservoir.R
@@ -68,10 +68,17 @@ RunModel_Reservoir <- function(InputsModel, RunOptions, Param) {
iPerTot <- seq(length(IndPerTot))
Vsim <- rep(0, length(IndPerTot))
Qsim_m3 <- Vsim
+ if (InputsModel$hasDiversion) {
+ Qdiv_m3 <- Vsim
+ }
# Time series volume and release calculation
for(i in iPerTot) {
Vsim[i] <- V0 + Qinflows_m3[i]
+ if (InputsModel$hasDiversion) {
+ Qdiv_m3[i] <- min(Vsim[i] + InputsModel$Qmin[IndPerTot[i]], InputsModel$Qdiv[IndPerTot[i]])
+ Vsim[i] <- Vsim[i] - Qdiv_m3[i]
+ }
Qsim_m3[i] <- min(Vsim[i], InputsModel$Qrelease[IndPerTot[i]])
Vsim[i] <- Vsim[i] - Qsim_m3[i]
if (Vsim[i] > Vmax) {
@@ -86,10 +93,16 @@ RunModel_Reservoir <- function(InputsModel, RunOptions, Param) {
iWarmUp <- seq(length(RunOptions$IndPeriod_WarmUp))
OutputsModel$RunOptions$WarmUpQsim_m3 <- Qsim_m3[iWarmUp]
OutputsModel$RunOptions$WarmUpVsim <- Vsim[iWarmUp]
+ if (InputsModel$hasDiversion) {
+ OutputsModel$RunOptions$WarmUpQdiv_m3 <- Qdiv_m3[iWarmUp]
+ }
}
iRun <- length(IndPerWarmUp) + seq(length(RunOptions$IndPeriod_Run))
OutputsModel$Qsim_m3 <- Qsim_m3[iRun]
OutputsModel$Vsim <- Vsim[iRun]
+ if (InputsModel$hasDiversion) {
+ OutputsModel$Qdiv_m3 <- Qdiv_m3[iRun]
+ }
OutputsModel$StateEnd$Reservoir <- list(V = Vsim[length(Vsim)])
class(OutputsModel) <- c("OutputsModelReservoir", class(OutputsModel))
return(OutputsModel)
diff --git a/tests/testthat/helper_1_RunModel.R b/tests/testthat/helper_1_RunModel.R
index 67033c4..ef1d0b7 100644
--- a/tests/testthat/helper_1_RunModel.R
+++ b/tests/testthat/helper_1_RunModel.R
@@ -51,6 +51,49 @@ setupRunModel <-
TempMean <- NULL
}
+ # Calibration parameters
+ ParamMichel <- list(
+ `54057` = c(
+ 0.781180650559296,
+ 74.4247133147015,
+ -1.26590474908235,
+ 0.96012365697022,
+ 2.51101785373787
+ ),
+ `54032` = c(
+ 0.992743594649893,
+ 1327.19309205366,
+ -0.505902143697436,
+ 7.91553615210291,
+ 2.03604525989572
+ ),
+ `54001` = c(
+ 1.03,
+ 24.7790862245877,
+ -1.90430150145153,
+ 21.7584023961971,
+ 1.37837837837838
+ ),
+ `54095` = c(
+ 256.844150254651,
+ 0.0650458497009288,
+ 57.523675209819,
+ 2.71809513102128
+ ),
+ `54002` = c(
+ 419.437754485522,
+ 0.12473266292168,
+ 13.0379482833606,
+ 2.12230907892238
+ ),
+ `54029` = c(
+ 219.203385553954,
+ 0.389211590110934,
+ 48.4242150713452,
+ 2.00300300300301
+ )
+ )
+
# set up inputs
if (!runInputsModel)
return(environment())
@@ -85,49 +128,5 @@ setupRunOptions <- function(InputsModel) {
RunOptions <- CreateRunOptions(InputsModel,
IndPeriod_WarmUp = IndPeriod_WarmUp,
IndPeriod_Run = IndPeriod_Run)
-
- # Calibration parameters
- ParamMichel <- list(
- `54057` = c(
- 0.781180650559296,
- 74.4247133147015,
- -1.26590474908235,
- 0.96012365697022,
- 2.51101785373787
- ),
- `54032` = c(
- 0.992743594649893,
- 1327.19309205366,
- -0.505902143697436,
- 7.91553615210291,
- 2.03604525989572
- ),
- `54001` = c(
- 1.03,
- 24.7790862245877,
- -1.90430150145153,
- 21.7584023961971,
- 1.37837837837838
- ),
- `54095` = c(
- 256.844150254651,
- 0.0650458497009288,
- 57.523675209819,
- 2.71809513102128
- ),
- `54002` = c(
- 419.437754485522,
- 0.12473266292168,
- 13.0379482833606,
- 2.12230907892238
- ),
- `54029` = c(
- 219.203385553954,
- 0.389211590110934,
- 48.4242150713452,
- 2.00300300300301
- )
- )
-
return(environment())
}
diff --git a/tests/testthat/test-RunModel_Reservoir.R b/tests/testthat/test-RunModel_Reservoir.R
index 76ac4ec..7417935 100644
--- a/tests/testthat/test-RunModel_Reservoir.R
+++ b/tests/testthat/test-RunModel_Reservoir.R
@@ -102,6 +102,15 @@ test_that("Calibration with ungauged node and reservoir filled by a diversion wo
test_that("Diversion on a reservoir works #146", {
e <- setupRunModel(runInputsModel = FALSE)
for (x in ls(e)) assign(x, get(x, e))
+ Qrelease <- data.frame(Dam = rep(3508465, length(DatesR)))
+ Param <- c(ParamMichel[names(ParamMichel) %in% griwrm$id], list(Dam = c(10E6, 1)))
+ e <- setupRunModel(runRunModel = FALSE,
+ griwrm = CreateGRiwrm(n_rsrvr),
+ Qrelease = Qrelease)
+ for (x in ls(e)) assign(x, get(x, e))
+ OM_resOnly <- RunModel(InputsModel,
+ RunOptions = RunOptions,
+ Param = Param)
nodes <- rbind(
n_rsrvr,
data.frame(
@@ -112,17 +121,17 @@ test_that("Diversion on a reservoir works #146", {
model = "Diversion"
)
)
- Qrelease <- data.frame(Dam = rep(3508465, length(DatesR)))
- Qobs2 <- Qrelease / 10
+ Qobs2 <- Qrelease * 0.1
g <- CreateGRiwrm(nodes)
e <- setupRunModel(griwrm = g,
runRunModel = FALSE,
Qobs2 = Qobs2,
Qrelease = Qrelease)
for (x in ls(e)) assign(x, get(x, e))
- Param <- c(ParamMichel[names(ParamMichel) %in% griwrm$id], list(Dam = c(10E6, 1)))
+
OM <- RunModel(InputsModel,
RunOptions = RunOptions,
Param = Param)
expect_true(max(OM$Dam$Vsim) - min(OM$Dam$Vsim) > 0)
+ expect_false(all(OM$Dam$Vsim == OM_resOnly$Dam$Vsim))
})
--
GitLab