diff --git a/dev/articles/V04_Closed-loop_regulated_withdrawal.html b/dev/articles/V04_Closed-loop_regulated_withdrawal.html index 20e4cad..4afcd34 100644 --- a/dev/articles/V04_Closed-loop_regulated_withdrawal.html +++ b/dev/articles/V04_Closed-loop_regulated_withdrawal.html @@ -582,7 +582,7 @@

The control logic function
 str(as.list(environment(fIrrigation)))
 #> List of 4
-#>  $ supervisor          :Classes 'Supervisor', 'environment' <environment: 0x56075aef4f00> 
+#>  $ supervisor          :Classes 'Supervisor', 'environment' <environment: 0x55a604cf81f8> 
 #>  $ irrigationObjective :'data.frame':    12 obs. of  4 variables:
 #>   ..$ month: num [1:12] 1 2 3 4 5 6 7 8 9 10 ...
 #>   ..$ 54001: num [1:12] 0 0 0.4 0.8 1 1.2 1.2 1 0.6 0 ...
diff --git a/dev/pkgdown.yml b/dev/pkgdown.yml
index 63a8cdd..0a6ac8b 100644
--- a/dev/pkgdown.yml
+++ b/dev/pkgdown.yml
@@ -16,7 +16,7 @@ articles:
   seinebasin/V05b_Open-loop_influenced_flow_calibration_GR6J: seinebasin/V05b_Open-loop_influenced_flow_calibration_GR6J.html
   V06_Modelling_regulated_diversion: V06_Modelling_regulated_diversion.html
   seinebasin/V06_Naturalised_flow_simulation: seinebasin/V06_Naturalised_flow_simulation.html
-last_built: 2024-09-23T10:27Z
+last_built: 2024-09-23T10:41Z
 urls:
   reference: https://inrae.github.io/airGRiwrm/reference
   article: https://inrae.github.io/airGRiwrm/articles
diff --git a/dev/reference/CreateInputsModel.GRiwrm.html b/dev/reference/CreateInputsModel.GRiwrm.html
index 0604071..705b09c 100644
--- a/dev/reference/CreateInputsModel.GRiwrm.html
+++ b/dev/reference/CreateInputsModel.GRiwrm.html
@@ -152,6 +152,7 @@ 
Creation of an InputsModel object for a airGRiwrm network  HypsoData = NULL,
   NLayers = 5,
   IsHyst = FALSE,
+  FUN_REGUL = NULL,
   ...
 )
@@ -249,6 +250,10 @@

Arguments

CemaNeige is used. See details of airGR::CreateRunOptions.

+
FUN_REGUL
+

List of functions for local regulation (See details)

+ +
...

used for compatibility with S3 methods

@@ -268,16 +273,36 @@

Details

requires named vector with the id of the sub-catchment as name item. If an argument is optional, only the column or the named item has to be provided.

See airGR::CreateInputsModel documentation for details concerning each input.

-

Number of rows of Precip, PotEvap, Qinf, Qmin, TempMean, TempMin, -TempMax must be the same of the length of DatesR (each row corresponds to -a time step defined in DatesR).

+

Number of rows of Precip, PotEvap, Qinf, Qmin, Qrelease, TempMean, +TempMin, TempMax must be the same of the length of DatesR (each row +corresponds to a time step defined in DatesR).

For examples of use see topics RunModel.GRiwrmInputsModel, RunModel_Reservoir, and RunModel.Supervisor.

For example of use of Direct Injection nodes, see vignettes "V03_Open-loop_influenced_flow" and "V04_Closed-loop_regulated_withdrawal".

For example of use of Diversion nodes, see example in RunModel.GRiwrmInputsModel topic and vignette -"V06_Modelling_regulated_diversion".

+"V06_Modelling_regulated_diversion".

+

The FUN_REGUL parameter

+ + +

FUN_REGUL argument is a named list of function that modify the node +InputsModel before sending it to the node's model. +This feature is useful for modifying data such as InputsModel$Qdiv or +InputsModel$Qrelease giving simulated flows already available from upstream +nodes. +Each item of the list has a name corresponding to the node on which the +function is applied. Each function must follow this interface: +function(InputsModel, RunOptions, OutputsModel, env) where the arguments are:

The functions embedded in FUN_REGUL should all return the argument +InputsModel after calculation.

+
+

See also

diff --git a/dev/reference/RunModel.Supervisor.html b/dev/reference/RunModel.Supervisor.html index 5cea9d4..813be28 100644 --- a/dev/reference/RunModel.Supervisor.html +++ b/dev/reference/RunModel.Supervisor.html @@ -212,14 +212,12 @@

Examples

# This release will be modified by the Supervisor # We initiate it with the natural flow for having a good initialization of the # model at the first time step of the running period -Qinf <- data.frame( +Qrelease <- data.frame( Dam = BasinsObs$`54095`$discharge_spec * griwrm$area[griwrm$id == "54095"] * 1E3 ) # InputsModel object -IM_severn <- CreateInputsModel(griwrm, DatesR, Precip, PotEvap, Qinf) -#> Warning: Use of the `Qinf` parameter for reservoir releases is deprecated, please use `Qrelease` instead. -#> Processing `Qrelease <- cbind(Qrelease, Qinf[, c("Dam"])`... +IM_severn <- CreateInputsModel(griwrm, DatesR, Precip, PotEvap, Qrelease = Qrelease) #> CreateInputsModel.GRiwrm: Processing sub-basin 54095... #> CreateInputsModel.GRiwrm: Processing sub-basin 54002... #> CreateInputsModel.GRiwrm: Processing sub-basin 54029... diff --git a/dev/reference/RunModel_Reservoir-5.png b/dev/reference/RunModel_Reservoir-5.png new file mode 100644 index 0000000..aac1ec4 Binary files /dev/null and b/dev/reference/RunModel_Reservoir-5.png differ diff --git a/dev/reference/RunModel_Reservoir-6.png b/dev/reference/RunModel_Reservoir-6.png new file mode 100644 index 0000000..505365a Binary files /dev/null and b/dev/reference/RunModel_Reservoir-6.png differ diff --git a/dev/reference/RunModel_Reservoir.html b/dev/reference/RunModel_Reservoir.html index aa374e0..c1c2170 100644 --- a/dev/reference/RunModel_Reservoir.html +++ b/dev/reference/RunModel_Reservoir.html @@ -230,7 +230,8 @@

Examples

# We propose to compute the constant released flow from # the median of the natural flow # The value is in m3 by time step (day) -Qrelease <- median(BasinObs$Qls, na.rm = TRUE) / 1000 * 86400 +(Qrelease <- median(BasinObs$Qls, na.rm = TRUE) / 1000 * 86400) +#> [1] 349920 # Formatting of reservoir released flow inputs for airGRiwrm (matrix or data.frame # with one column by node and node IDs as column names) @@ -239,9 +240,7 @@

Examples

InputsModel <- CreateInputsModel(griwrm, DatesR = BasinObs$DatesR, Precip = Precip, PotEvap = PotEvap, - Qinf = Qrelease) -#> Warning: Use of the `Qinf` parameter for reservoir releases is deprecated, please use `Qrelease` instead. -#> Processing `Qrelease <- cbind(Qrelease, Qinf[, c("Reservoir"])`... + Qrelease = Qrelease) #> CreateInputsModel.GRiwrm: Processing sub-basin L0123001... #> CreateInputsModel.GRiwrm: Processing sub-basin Reservoir... @@ -267,15 +266,11 @@

Examples

RunOptions = RunOptions, Obs = Qobs) -# preparation of CalibOptions object -CalibOptions <- CreateCalibOptions(InputsModel) -#> Warning: The node 'Reservoir' which uses `RunModel_Reservoir` must have its parameters fixed: -#> You can either fix these parameters afterward by using the command: -#> `CalibOptions[['Reservoir']]$FixedParam <- c(Vmax, celerity)` -#> Or by calling `CreateCalibOptions(InputsModel, FixedParam = list('Reservoir' = c(Vmax, celerity)))` - -# Parameters of RunModel_Reservoir must be fixed -CalibOptions[["Reservoir"]]$FixedParam <- c(Vmax = 30E6, celerity = 0.5) +# preparation of CalibOptions object with fixed parameters for the reservoir +Vmax <- 30E6 +CalibOptions <- + CreateCalibOptions(InputsModel, + FixedParam = list(Reservoir = c(Vmax = 30E6, celerity = 0.5))) OC <- Calibration( InputsModel = InputsModel, @@ -326,6 +321,92 @@

Examples

# The plot for the reservoir can also be plotted alone plot(OutputsModel$Reservoir, Qobs = Qobs[, "Reservoir"]) + +####################################################### +# Daily time step simulation of a reservoir tracking # +# an objective filling curve using a local regulation # +####################################################### + +# The objective here is to simulate the same reservoir as above +# but with new rules: +# - A minimum flow downstream the reservoir defined as: +(Qmin <- Qrelease[1,] / 2) +#> [1] 174960 +# - A maximum release flow due to reservoir outlet limitation +(Qmax <- Qrelease[1,] * 5) +#> [1] 1749600 +# - An annual objective filling curve managing floods and droughts by +# trying to keep the reservoir volume between 10 and 20 Mm3: +Vobj <- approx(c(1, 120, 300, 366), + c(20E6, 20E6, 10E6, 20E6), + seq(366)) +plot(Vobj, type = "l", col = "red", lty = 2) + + +# The regulation function takes InputsModel of the reservoir node and the +# global GRiwrm OutputsModel as arguments and returns a modified +# InputsModel used by RunModel_Reservoir afterward +fun_factory_Regulation_Reservoir <- function(Vini, Vobj, Qmin, Qmax, Vmax) { + function(InputsModel, RunOptions, OutputsModel, env) { + # Release flow time series initialisation + Qrelease <- rep(0, length(InputsModel$DatesR)) + # Build inflows time series from upstream Qsim (warmup & run) + Qinflows <- Qrelease + IPR_all <- c(RunOptions$IndPeriod_WarmUp, RunOptions$IndPeriod_Run) + Qinflows[IPR_all] <- c(OutputsModel$L0123001$RunOptions$WarmUpQsim_m3, + OutputsModel$L0123001$Qsim_m3) + # Reservoir volume initialisation + V <- Vini + # Loop over simulation time steps (warmup & run periods) + for(ts in IPR_all) { + # Update reservoir volume with inflows + V <- V + Qinflows[ts] + # Rule #1: follow the objective filling curve (lower priority) + j <- as.numeric(format(InputsModel$DatesR[ts], "%j")) + Vobj_ts <- approx(Vobj, xout = j)$y + Qrelease[ts] <- V - Vobj_ts + # Rule #2: Release cannot be less than Qmin + Qrelease[ts] <- max(Qmin, Qrelease[ts]) + # Rule #3: Release cannot be more than Qmax + Qrelease[ts] <- min(Qmax, Qrelease[ts]) + # Update reservoir volume after release + V <- V - Qrelease[ts] + # Rule #4: hard constraints on the reservoir (full or empty?) + if (V < 0) { + Qrelease[ts] <- Qrelease[ts] + V + V <- 0 + } + V <- min(V, Vmax) + } + InputsModel$Qrelease <- Qrelease + return(InputsModel) + } +} +# A call to fun_factory_Regulation_Reservoir returns the regulation +# function with the parameters Qmin, Qmax, Vobj enclosed in the environment +# of the function +Regulation_Reservoir <- + fun_factory_Regulation_Reservoir(RunOptions$Reservoir$IniStates, Vobj, Qmin, Qmax, Vmax) + +# Then we need to update InputsModel in order to take into account the regulation +# function instead of predefined Qrelease in the previous study case +IM_reg <- CreateInputsModel(griwrm, + DatesR = BasinObs$DatesR, + Precip = Precip, + PotEvap = PotEvap, + Qrelease = Qrelease, + FUN_REGUL = list(Reservoir = Regulation_Reservoir)) +#> CreateInputsModel.GRiwrm: Processing sub-basin L0123001... +#> CreateInputsModel.GRiwrm: Processing sub-basin Reservoir... + +# And we can finally run the simulation! +OM_reg <- RunModel(IM_reg, RunOptions, Param) +#> RunModel.GRiwrmInputsModel: Processing sub-basin L0123001... +#> RunModel.GRiwrmInputsModel: Processing sub-basin Reservoir... + +# And plot the new result +plot(OM_reg$Reservoir) +
diff --git a/dev/reference/mermaid.html b/dev/reference/mermaid.html index ec816c6..73ae23c 100644 --- a/dev/reference/mermaid.html +++ b/dev/reference/mermaid.html @@ -226,11 +226,11 @@

Examples

#> [1] "https://mermaid.ink/img/pako:eNqrVkrOT0lVslJKy8kvT85ILCpR8AmKyVNQcFTQ1bVTcFLSUcpNLcpNzExRsqpWKslIzQUpTklNSyzNKVGqrQUAjIcUfg?type=png" f <- mermaid(diagram) f -#> [1] "/tmp/Rtmp8s6lXi/a0a48ce294b1a6f989506ce145809601.png" +#> [1] "/tmp/RtmpyQ3KUW/a0a48ce294b1a6f989506ce145809601.png" # For displaying the diagram in Rmarkdown document knitr::include_graphics(mermaid(diagram)) -#> [1] "/tmp/Rtmp8s6lXi/a0a48ce294b1a6f989506ce145809601.png" +#> [1] "/tmp/RtmpyQ3KUW/a0a48ce294b1a6f989506ce145809601.png" #> attr(,"class") #> [1] "knit_image_paths" "knit_asis" diff --git a/dev/reference/plot.OutputsModelReservoir-5.png b/dev/reference/plot.OutputsModelReservoir-5.png new file mode 100644 index 0000000..aac1ec4 Binary files /dev/null and b/dev/reference/plot.OutputsModelReservoir-5.png differ diff --git a/dev/reference/plot.OutputsModelReservoir-6.png b/dev/reference/plot.OutputsModelReservoir-6.png new file mode 100644 index 0000000..505365a Binary files /dev/null and b/dev/reference/plot.OutputsModelReservoir-6.png differ diff --git a/dev/reference/plot.OutputsModelReservoir.html b/dev/reference/plot.OutputsModelReservoir.html index 3ac0006..f9b6e47 100644 --- a/dev/reference/plot.OutputsModelReservoir.html +++ b/dev/reference/plot.OutputsModelReservoir.html @@ -193,7 +193,8 @@

Examples

# We propose to compute the constant released flow from # the median of the natural flow # The value is in m3 by time step (day) -Qrelease <- median(BasinObs$Qls, na.rm = TRUE) / 1000 * 86400 +(Qrelease <- median(BasinObs$Qls, na.rm = TRUE) / 1000 * 86400) +#> [1] 349920 # Formatting of reservoir released flow inputs for airGRiwrm (matrix or data.frame # with one column by node and node IDs as column names) @@ -202,9 +203,7 @@

Examples

InputsModel <- CreateInputsModel(griwrm, DatesR = BasinObs$DatesR, Precip = Precip, PotEvap = PotEvap, - Qinf = Qrelease) -#> Warning: Use of the `Qinf` parameter for reservoir releases is deprecated, please use `Qrelease` instead. -#> Processing `Qrelease <- cbind(Qrelease, Qinf[, c("Reservoir"])`... + Qrelease = Qrelease) #> CreateInputsModel.GRiwrm: Processing sub-basin L0123001... #> CreateInputsModel.GRiwrm: Processing sub-basin Reservoir... @@ -230,15 +229,11 @@

Examples

RunOptions = RunOptions, Obs = Qobs) -# preparation of CalibOptions object -CalibOptions <- CreateCalibOptions(InputsModel) -#> Warning: The node 'Reservoir' which uses `RunModel_Reservoir` must have its parameters fixed: -#> You can either fix these parameters afterward by using the command: -#> `CalibOptions[['Reservoir']]$FixedParam <- c(Vmax, celerity)` -#> Or by calling `CreateCalibOptions(InputsModel, FixedParam = list('Reservoir' = c(Vmax, celerity)))` - -# Parameters of RunModel_Reservoir must be fixed -CalibOptions[["Reservoir"]]$FixedParam <- c(Vmax = 30E6, celerity = 0.5) +# preparation of CalibOptions object with fixed parameters for the reservoir +Vmax <- 30E6 +CalibOptions <- + CreateCalibOptions(InputsModel, + FixedParam = list(Reservoir = c(Vmax = 30E6, celerity = 0.5))) OC <- Calibration( InputsModel = InputsModel, @@ -289,6 +284,92 @@

Examples

# The plot for the reservoir can also be plotted alone plot(OutputsModel$Reservoir, Qobs = Qobs[, "Reservoir"]) + +####################################################### +# Daily time step simulation of a reservoir tracking # +# an objective filling curve using a local regulation # +####################################################### + +# The objective here is to simulate the same reservoir as above +# but with new rules: +# - A minimum flow downstream the reservoir defined as: +(Qmin <- Qrelease[1,] / 2) +#> [1] 174960 +# - A maximum release flow due to reservoir outlet limitation +(Qmax <- Qrelease[1,] * 5) +#> [1] 1749600 +# - An annual objective filling curve managing floods and droughts by +# trying to keep the reservoir volume between 10 and 20 Mm3: +Vobj <- approx(c(1, 120, 300, 366), + c(20E6, 20E6, 10E6, 20E6), + seq(366)) +plot(Vobj, type = "l", col = "red", lty = 2) + + +# The regulation function takes InputsModel of the reservoir node and the +# global GRiwrm OutputsModel as arguments and returns a modified +# InputsModel used by RunModel_Reservoir afterward +fun_factory_Regulation_Reservoir <- function(Vini, Vobj, Qmin, Qmax, Vmax) { + function(InputsModel, RunOptions, OutputsModel, env) { + # Release flow time series initialisation + Qrelease <- rep(0, length(InputsModel$DatesR)) + # Build inflows time series from upstream Qsim (warmup & run) + Qinflows <- Qrelease + IPR_all <- c(RunOptions$IndPeriod_WarmUp, RunOptions$IndPeriod_Run) + Qinflows[IPR_all] <- c(OutputsModel$L0123001$RunOptions$WarmUpQsim_m3, + OutputsModel$L0123001$Qsim_m3) + # Reservoir volume initialisation + V <- Vini + # Loop over simulation time steps (warmup & run periods) + for(ts in IPR_all) { + # Update reservoir volume with inflows + V <- V + Qinflows[ts] + # Rule #1: follow the objective filling curve (lower priority) + j <- as.numeric(format(InputsModel$DatesR[ts], "%j")) + Vobj_ts <- approx(Vobj, xout = j)$y + Qrelease[ts] <- V - Vobj_ts + # Rule #2: Release cannot be less than Qmin + Qrelease[ts] <- max(Qmin, Qrelease[ts]) + # Rule #3: Release cannot be more than Qmax + Qrelease[ts] <- min(Qmax, Qrelease[ts]) + # Update reservoir volume after release + V <- V - Qrelease[ts] + # Rule #4: hard constraints on the reservoir (full or empty?) + if (V < 0) { + Qrelease[ts] <- Qrelease[ts] + V + V <- 0 + } + V <- min(V, Vmax) + } + InputsModel$Qrelease <- Qrelease + return(InputsModel) + } +} +# A call to fun_factory_Regulation_Reservoir returns the regulation +# function with the parameters Qmin, Qmax, Vobj enclosed in the environment +# of the function +Regulation_Reservoir <- + fun_factory_Regulation_Reservoir(RunOptions$Reservoir$IniStates, Vobj, Qmin, Qmax, Vmax) + +# Then we need to update InputsModel in order to take into account the regulation +# function instead of predefined Qrelease in the previous study case +IM_reg <- CreateInputsModel(griwrm, + DatesR = BasinObs$DatesR, + Precip = Precip, + PotEvap = PotEvap, + Qrelease = Qrelease, + FUN_REGUL = list(Reservoir = Regulation_Reservoir)) +#> CreateInputsModel.GRiwrm: Processing sub-basin L0123001... +#> CreateInputsModel.GRiwrm: Processing sub-basin Reservoir... + +# And we can finally run the simulation! +OM_reg <- RunModel(IM_reg, RunOptions, Param) +#> RunModel.GRiwrmInputsModel: Processing sub-basin L0123001... +#> RunModel.GRiwrmInputsModel: Processing sub-basin Reservoir... + +# And plot the new result +plot(OM_reg$Reservoir) +