diff --git a/R/ogs6_ensemble.R b/R/ogs6_ensemble.R
index 09eae5a8700049050a897535166245ebe827724d..76ea85952d30fe4e8bdf99b17901e743fd8e7bac 100644
--- a/R/ogs6_ensemble.R
+++ b/R/ogs6_ensemble.R
@@ -286,10 +286,9 @@ OGS6_Ensemble <- R6::R6Class(
val <- eval(parse(text = self$dp_parameters[[i]]))
val <- as.numeric(val)
- val_vec <- vapply(private$.parameter_percs[[i]], function(x){
+ val_vec <- lapply(private$.parameter_percs[[i]], function(x){
val + (val * (x / 100))
- },
- FUN.VALUE = numeric(length(val)))
+ })
private$.parameter_values <- c(self$parameter_values,
list(val_vec))
diff --git a/r2ogs6.Rproj b/r2ogs6.Rproj
index 4ba229d3862d841ed6e19855704fa9287bef7626..7cca3129cdb7302da383fed22c8a428cac9b3356 100644
--- a/r2ogs6.Rproj
+++ b/r2ogs6.Rproj
@@ -18,4 +18,5 @@ BuildType: Package
PackageUseDevtools: Yes
PackageInstallArgs: --no-multiarch --with-keep.source
PackageBuildArgs: --no-build-vignettes
+PackageCheckArgs: --no-build-vignettes
PackageRoxygenize: rd,collate,namespace
diff --git a/vignettes/ensemble_workflow_vignette.Rmd b/vignettes/ensemble_workflow_vignette.Rmd
index 4cacf5bd75db36747ecfc6880ea24ce58bbebaed..5b598b24e9cf0d6104037e1aec3e0076eca255b1 100644
--- a/vignettes/ensemble_workflow_vignette.Rmd
+++ b/vignettes/ensemble_workflow_vignette.Rmd
@@ -9,6 +9,7 @@ vignette: >
```{r, include = FALSE}
knitr::opts_chunk$set(
+ eval = nzchar(Sys.getenv("r2ogs6_ensemble_guide_eval")),
fig.width = 6,
fig.align = "center",
collapse = TRUE,
@@ -55,7 +56,7 @@ First, we create a simulation object to base our ensemble on and read in the `.p
```{r}
# Change this to fit your system
-testdir_path <- system.file("extdata/test_tempdirs/", package = "r2ogs6")
+testdir_path <- system.file("extdata/test_tempdirs", package = "r2ogs6")
sim_path <- paste0(testdir_path, "/axisym_theis_sim")
ogs6_obj <- OGS6$new(sim_name = "axisym_theis",
@@ -223,11 +224,11 @@ We will consider the following parameters for our sensitivity analysis:
### Setup
-First, we create a simulation object to base our ensemble on and read in the `.prj` file.
+First, we create a simulation object to base our ensemble on and read in the `.prj` file. This time we want to specify that an output file only gets written at the last timestep.
```{r}
# Change this to fit your system
-testdir_path <- system.file("extdata/test_tempdirs/", package = "r2ogs6")
+testdir_path <- system.file("extdata/test_tempdirs", package = "r2ogs6")
sim_path <- paste0(testdir_path, "/theis_sim")
ogs6_obj <- OGS6$new(sim_name = "theis",
@@ -238,6 +239,9 @@ prj_path <- system.file("examples/Theis_well_pumping/benchmark_files/",
"theis.prj", package = "r2ogs6")
read_in_prj(ogs6_obj, prj_path)
+
+# Increase each_steps
+ogs6_obj$time_loop$output$timesteps$pair$each_steps <- 200
```
### Multi-parameter sequential run
@@ -251,13 +255,13 @@ ogs6_ens_theis_2 <-
ogs6_obj = ogs6_obj,
parameters =
list(
- per = list(
+ per = list(ogs6_obj$parameters[[3]]$values,
+ values = percentages),
+ por = list(ogs6_obj$parameters[[2]]$value,
+ values = percentages),
+ slo = list(
ogs6_obj$media[[1]]$phases[[1]]$properties[[1]]$independent_variable[[2]]$slope,
- values = percentages),
- por = list(ogs6_obj$media[[1]]$properties[[2]]$value,
- c(-50, -10, -1, 0)),
- slo = list(ogs6_obj$media[[1]]$properties[[4]]$value,
- values = percentages)
+ values = percentages)
),
sequential_mode = TRUE
)
@@ -266,12 +270,29 @@ exit_codes <- ogs6_ens_theis_2$run_simulation()
paste(exit_codes)
```
-When the simulations have run, we can extract and plot the results just like before.
+When the simulations have run, we can extract and plot the results like before. To avoid cluttering the plot, we only extract the `pressure` values for a single line. For this, we get the IDs of all points on the x axis.
```{r}
+# Extract point ids
+get_point_ids_x <- function(points){
+ x_axis_ids <- numeric()
+
+ for(i in seq_len(dim(points)[[1]])) {
+ if (points[i, ][[2]] == 0 && points[i, ][[3]] == 0) {
+ x_axis_ids <- c(x_axis_ids, (i - 1))
+ }
+ }
+
+ return(x_axis_ids)
+}
+
+point_ids_x <- get_point_ids_x(
+ ogs6_ens_theis_2$ensemble[[1]]$pvds[[1]]$OGS6_vtus[[1]]$points)
+
# Get combined dataframe
per_por_slo_df <-
ogs6_ens_theis_2$get_point_data(
+ point_ids = point_ids_x,
keys = c("pressure"),
start_at_timestep = ogs6_ens_theis_2$ensemble[[1]]$pvds[[1]]$last_timestep)
```
@@ -288,6 +309,82 @@ ggplot(per_por_slo_df,
facet_grid(cols = vars(name))
```
+
+Let's take a closer look at `permeability`.
+
+```{r}
+per_df <- subset(per_por_slo_df, name == "per")
+
+# Make plot
+ggplot(per_df,
+ aes(x = x,
+ y = pressure)) +
+ geom_point(aes(color = as.factor(perc))) +
+ xlab("x point coordinate") +
+ labs(color = "%")
+```
+
+Maybe we want to try and use a logarithmic approach for `slope`. This won't work with the built-in functionality of `OGS6_Ensemble` so we'll set up our Ensemble a little differently.
+
+```{r}
+# Calculate log value
+log_val <- log(as.numeric(
+ ogs6_obj$media[[1]]$phases[[1]]$properties[[1]]$independent_variable[[2]]$slope),
+ base = 10)
+
+# Apply changes to log value
+log_vals <- vapply(percentages, function(x){
+ log_val + (log_val * (x / 100))
+}, FUN.VALUE = numeric(1))
+
+# Transfer back to non-logarithmic scale
+back_transf_vals <- 10^log_vals
+
+# Set up new ensemble
+ogs6_ens_slo <-
+ OGS6_Ensemble$new(
+ ogs6_obj = ogs6_obj,
+ parameters =
+ list(
+ slo = list(
+ ogs6_obj$media[[1]]$phases[[1]]$properties[[1]]$independent_variable[[2]]$slope,
+ values = back_transf_vals)
+ ),
+ percentages_mode = FALSE,
+ sequential_mode = TRUE
+ )
+
+exit_codes <- ogs6_ens_slo$run_simulation()
+```
+
+Let's check if we can observe any influence of `slope` on `pressure` now.
+
+```{r}
+# Filter point ids
+point_ids_x <- get_point_ids_x(
+ ogs6_ens_slo$ensemble[[1]]$pvds[[1]]$OGS6_vtus[[1]]$points)
+
+# Get combined dataframe
+slo_df <-
+ ogs6_ens_slo$get_point_data(
+ point_ids = point_ids_x,
+ keys = c("pressure"),
+ start_at_timestep = ogs6_ens_slo$ensemble[[1]]$pvds[[1]]$last_timestep)
+
+# Supply percentages manually since we couldn't use `percentages_mode`
+percs <- vapply(slo_df$sim_id,
+ function(x){percentages[[x]]},
+ FUN.VALUE = numeric(1))
+
+ggplot(slo_df,
+ aes(x = x,
+ y = pressure)) +
+ geom_point(aes(color = as.factor(percs))) +
+ xlab("x point coordinate") +
+ labs(color = "%")
+```
+
+
## Summary
The `OGS6_Ensemble` class is a useful tool to set up ensemble runs for sensitivity analyses. In this vignette, we learned how to create `OGS6_Ensemble` objects. We looked at how the parameters `sequential_mode` and `percentages_mode` influence how our ensemble object is initialised. We started simulations via `OGS6_Ensemble$run_simulation()` and extracted information from the output files to plot them.