Skip to content
GitLab
Select Git revision
  • 181-simulation-water-losses-on-routed-flows-based-on-gr-groundwater-exchange-function
  • dev default
  • 183-createcaliboptions-crash-when-all-parameters-are-fixed
  • 180-write-a-vignette-with-exemple-of-use-of-runmodel-griwrmoutputsmodel
  • master
  • v0.7.0
  • v0.6.2
  • v0.6.1
  • v0.5.0
  • v0.4.0
  • v0.3.0
Find file BlameHistoryPermalink
RunModel_Lag.R 3.97 KiB
1 
RunModel_Lag <- function(InputsModel, RunOptions, Param) {
2 
  NParam <- 1
3
4 
  ##Arguments_check
5 
  if (!inherits(InputsModel, "InputsModel")) {
6 
    stop("'InputsModel' must be of class 'InputsModel'")
7 
  }
8 
  if (!inherits(InputsModel, "SD")) {
9 
    stop("'InputsModel' must be of class 'SD'")
10 
  }
11 
  if (!inherits(RunOptions, "RunOptions")) {
12 
    stop("'RunOptions' must be of class 'RunOptions'")
13 
  }
14 
  if (!is.vector(Param) | !is.numeric(Param)) {
15 
    stop("'Param' must be a numeric vector")
16 
  }
17 
  if (sum(!is.na(Param)) != NParam) {
18 
    stop(paste("'Param' must be a vector of length", NParam, "and contain no NA"))
19 
  }
20 
  if (is.null(InputsModel$OutputsModel)) {
21 
    stop(
22 
      "'InputsModel' should contain an 'OutputsModel' key containing the output of the runoff of the downstream subcatchment"
23 
    )
24 
  }
25 
  if (is.null(InputsModel$OutputsModel$Qsim)) {
26 
    stop(
27 
      "'InputsModel$OutputsModel' should contain a key 'Qsim' containing the output of the runoff of the downstream subcatchment"
28 
    )
29 
  }
30 
  if (sum(!is.na(InputsModel$OutputsModel$Qsim)) != length(RunOptions$IndPeriod_Run)) {
31 
    stop(
32 
      "'InputsModel$OutputsModel$Qim' should have the same lenght as 'RunOptions$IndPeriod_Run' and contain no NA"
33 
    )
34 
  }
35
36 
  OutputsModel <- InputsModel$OutputsModel
37 
  OutputsModel$QsimDown <- OutputsModel$Qsim
38
39 
  if (inherits(InputsModel, "hourly")) {
40 
    TimeStep <- 60 * 60
41 
  } else if (inherits(InputsModel, "daily")) {
42 
    TimeStep <- 60 * 60 * 24
43 
  } else {
44 
    stop("'InputsModel' should be of class \"daily\" or \"hourly\"")
45 
  }
46
47 
  # propagation time from upstream meshes to outlet
48 
  PT <- InputsModel$LengthHydro / Param[1L] / TimeStep
49 
  HUTRANS <- rbind(1 - (PT - floor(PT)), PT - floor(PT))
50
51 
  NbUpBasins <- length(InputsModel$LengthHydro)
52 
  LengthTs <- length(OutputsModel$QsimDown)
53 
  OutputsModel$Qsim <- OutputsModel$QsimDown * InputsModel$BasinAreas[length(InputsModel$BasinAreas)] * 1e3
54
55 
  IniSD <- RunOptions$IniStates[grep("SD", names(RunOptions$IniStates))]
56 
  if (length(IniSD) > 0) {
57 
    if (sum(floor(PT)) + NbUpBasins != length(IniSD)) {
58 
      stop(
59 
        sprintf(
60 
          "SD initial states has a length of %i and a length of %i is required",
61 
          length(IniSD),
62 
          sum(floor(PT)) + NbUpBasins
63 
        )
64 
      )
65 
    }
66 
    IniStates <- lapply(seq_len(NbUpBasins), function(x) {
67 
      iStart <- 1
68 
      if (x > 1) {
69 
        iStart <- iStart + sum(floor(PT[1:x - 1]) + 1)
70 
      }
7172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112
IniSD[iStart:(iStart + PT[x])] }) } else { IniStates <- lapply(seq_len(NbUpBasins), function(x) { rep(0, floor(PT[x] + 1)) }) } #message("Initstates: ",paste(IniStates, collapse = ", ")) for (upstream_basin in seq_len(NbUpBasins)) { Qupstream <- c(IniStates[[upstream_basin]], InputsModel$Qupstream[RunOptions$IndPeriod_Run, upstream_basin]) if (!is.na(InputsModel$BasinAreas[upstream_basin])) { # Upstream flow with area needs to be converted to m3 by time step Qupstream <- Qupstream * InputsModel$BasinAreas[upstream_basin] * 1e3 } #message("Qupstream[", upstream_basin, "]: ", paste(Qupstream, collapse = ", ")) OutputsModel$Qsim <- OutputsModel$Qsim + Qupstream[2:(1 + LengthTs)] * HUTRANS[1, upstream_basin] + Qupstream[1:LengthTs] * HUTRANS[2, upstream_basin] } # Warning for negative flows if (any(OutputsModel$Qsim < 0)) { warning(length(which(OutputsModel$Qsim < 0)), " time steps with negative flow, set to zero.") OutputsModel$Qsim[OutputsModel$Qsim < 0] <- 0 } # Convert back Qsim to mm OutputsModel$Qsim <- OutputsModel$Qsim / sum(InputsModel$BasinAreas, na.rm = TRUE) / 1e3 #message("Qsim: ",paste(OutputsModel$Qsim, collapse = ", ")) if ("StateEnd" %in% RunOptions$Outputs_Sim) { OutputsModel$StateEnd$SD <- lapply(seq(NbUpBasins), function(x) { lastTS <- RunOptions$IndPeriod_Run[length(RunOptions$IndPeriod_Run)] InputsModel$Qupstream[(lastTS - floor(PT[x])):lastTS, x] }) #message("StateEnd: ",paste(OutputsModel$StateEnd$SD, collapse = ", ")) } return(OutputsModel) }