From 68225767e2b257b2ecdcd92bb73c5b7ba1bac74e Mon Sep 17 00:00:00 2001
From: Dorchies David <david.dorchies@inrae.fr>
Date: Wed, 29 Jun 2022 09:25:14 +0200
Subject: [PATCH] fix(RunModel_LLR): wrong parameter transformation
Previous one was from Munier but the model is not the same as Munier one.
This formulation is from https://doi.org/10.5194/nhess-14-2899-2014
Refs HYCAR-Hydro/airgr#152, HYCAR-Hydro/airgr#153
---
R/RunModel_LLR.R | 21 ++++++++++-----------
1 file changed, 10 insertions(+), 11 deletions(-)
diff --git a/R/RunModel_LLR.R b/R/RunModel_LLR.R
index c11b12d..a37601b 100644
--- a/R/RunModel_LLR.R
+++ b/R/RunModel_LLR.R
@@ -68,19 +68,19 @@ RunModel_LLR <- function(InputsModel, RunOptions, Param, QcontribDown) {
TParam <- Param[1L]
KParam <- Param[2L]
- ## propagation time from upstream meshes to outlet
- PT <- ifelse(InputsModel$LengthHydro == 0, 0, InputsModel$LengthHydro/max(InputsModel$LengthHydro)) *
- (TParam + KParam) - ifelse(InputsModel$LengthHydro == 0, 0,
- sqrt(InputsModel$LengthHydro/max(InputsModel$LengthHydro))) * KParam
- PT[PT < 0] <- 0
- PK <- ifelse(InputsModel$LengthHydro == 0,
- 0,
- sqrt(InputsModel$LengthHydro/max(InputsModel$LengthHydro)) * KParam)
+ PT <- InputsModel$LengthHydro * 1e3 / Param[1L] / RunOptions$FeatFUN_MOD$TimeStep
+ PK <- KParam * PT
C0 <- ifelse(PK ==0, 0, exp(-1/PK))
C1 <- ((PK/1) * (1-C0))-C0
C2 <- 1 - (PK/1) * (1-C0)
+ ## Lag set up
+
+ Lag <- function(t, Qroute, Qupstream){
+ Q <- c(Qupstream[1 : floor(t)], Qroute[1 : (length(Qroute) - floor(t))])
+ return(Q)
+ }
## set up initial states
IniSD <- RunOptions$IniStates[grep("SD", names(RunOptions$IniStates))]
if (length(IniSD) > 0) {
@@ -135,10 +135,9 @@ RunModel_LLR <- function(InputsModel, RunOptions, Param, QcontribDown) {
Qroute_t1 <- Qroute[q_upstream - 1]
}
Qroute[q_upstream] <- C0[upstream_basin] * Qroute_t1 + C1[upstream_basin] * Qupstream_m3[q_upstream - 1] + C2[upstream_basin] * Qupstream_m3[q_upstream]
- Qllr[q_upstream + PT[upstream_basin]] <- Qroute[q_upstream]
}
- Qsim_m3[, upstream_basin] <- Qllr[(1 +length(Qllr) - length(IndPeriod1)) :
- (length(IndPeriod1) + (length(Qllr) - length(IndPeriod1)))]
+ Qllr <- Lag(PT[upstream_basin], Qroute, Qupstream_m3)[IndPeriod1]
+ Qsim_m3[, upstream_basin] <- Qllr
}
Qsim_m3 <- cbind(QsimDown_input, Qsim_m3)
--
GitLab