Dev

[mrustl]

No description provided.

[claude]

Potential race condition on the shared master file when parallel workers target the same site.

The master path is derived only from paths$dir_input and scenario_prefix, so both future_lapply workers of the same site (different storage_type) point to the same <site>/input/o_master.h5. First call in each worker (there are two, one per storage type) will:

1. fs::file_copy(base.h5 → o_master.h5, overwrite = TRUE)
2. hdf5r::H5File$new(mp, mode = "a") — opens the shared master in append mode

If both prepare_master calls interleave, worker B's file_copy(overwrite=TRUE) can stomp on the master file A currently has open in append mode, and/or the second H5File$new(..., mode="a") may fail due to HDF5 file locking. The vignettes in workflow_optimisation.Rmd fan out over expand.grid(site, type), so both workers of a site race here.

Suggested fixes (in order of decreasing intrusiveness):

• Include a per-worker discriminator in the master filename (e.g. sprintf("%s_master_%s.h5", scenario_prefix, Sys.getpid())), so each worker has its own; or
• Serialise creation via a lockfile (fs::file_touch + while(...) polling) — brittle;
• Document that make_swale_runner() is not safe to reuse a dir_input across concurrent workers.

The one-master-per-worker approach doubles disk use temporarily but eliminates the race entirely.

[claude]

Dead ceiling() on an already-integer column.

timeseries_et$time was computed above (line 34–36) as as.integer(difftime(..., units = "hours")), which is already an integer. ceiling() on integer values is a no-op — this line does nothing (compare with line 53–54 for rain, where time is a double and ceiling actually rounds up).

Two options:

• If the intent was to round up after as.integer's truncation, replace the earlier as.integer(...) with as.double(...) so the truncation doesn't happen first, then ceiling here matters.
• If the intent is a no-op (ET is daily and always integer hours by construction, so nothing to round), drop lines 40–41 entirely.

[claude]

evap_min_scenario becomes character(0) when every evap_all entry is NA, and the title collapses to character(0).

When no scenarios are valid AND every value in evap_all is NA (pathological but possible for a fully-failed batch), the fallback ref_idx <- seq_along(evap_all) fires, then which.min(evap_all[ref_idx]) returns integer(0) (min on all-NA), so ref_idx[integer(0)] = integer(0) and evap_min_scenario = character(0).

At line 283–291 the German/English branch does paste0(..., evap_min_scenario, ...). paste0 with a character(0) argument returns character(0), which then flows into title <- paste0(title, "\n", character(0)) → character(0). Passing a character(0) title to ggplot2::labs() produces an unhelpful downstream error rather than the intended graceful degradation.

Suggested guard (minimal patch, right at the definition):

Suggested change

	evap_all <- simulation_results_optimisation[["element.WB_Evapotranspiration_"]]
	valid_mask <- !is.na(simulation_results_optimisation$n_overflows) &
	simulation_results_optimisation$n_overflows <= x_int & !is.na(evap_all)
	ref_idx <- if (any(valid_mask)) which(valid_mask) else seq_along(evap_all)
	evap_min <- suppressWarnings(min(evap_all[ref_idx], na.rm = TRUE))
	evap_min_scenario <- simulation_results_optimisation$scenario_name[
	ref_idx[which.min(evap_all[ref_idx])]]
	# Share of scenarios meeting the validity criterion (n_overflows <= x),
	# appended to the auto-generated title (a plotly-safe place -- ggplotly
	# drops ggplot subtitles).
	valid_pct <- round(100 * mean(
	simulation_results_optimisation$n_overflows <= x_int, na.rm = TRUE))
	share_txt <- switch(lang,
	de = paste0(valid_pct, " % mit <= ", x_int, " \u00dcberl\u00e4ufen"),
	en = paste0(valid_pct, " % with <= ", x_int, " overflows"))
	if (is.null(title)) {
	title <- paste0(def_title, " (", share_txt, ")")
	if (y_var == "cost_per_evap_pct") {
	# name the reference of the marginal measure in the title
	title <- paste0(title, "\n", switch(lang,
	de = paste0("Referenz: minimale Evapotranspiration der g\u00fcltigen ",
	"Szenarien (<= ", x_int, " \u00dcberl\u00e4ufe): ",
	round(evap_min, 1), " % (Szenario ",
	evap_min_scenario, ")"),
	en = paste0("Reference: minimum evapotranspiration of the valid ",
	"scenarios (<= ", x_int, " overflows): ",
	round(evap_min, 1), " % (scenario ",
	evap_min_scenario, ")")))
	}
	}
	evap_all <- simulation_results_optimisation[["element.WB_Evapotranspiration_"]]
	valid_mask <- !is.na(simulation_results_optimisation$n_overflows) &
	simulation_results_optimisation$n_overflows <= x_int & !is.na(evap_all)
	ref_idx <- if (any(valid_mask)) which(valid_mask) else seq_along(evap_all)
	evap_min <- suppressWarnings(min(evap_all[ref_idx], na.rm = TRUE))
	wm <- which.min(evap_all[ref_idx])
	evap_min_scenario <- if (length(wm) == 0L) {
	NA_character_
	} else {
	simulation_results_optimisation$scenario_name[ref_idx[wm]]
	}
	# Share of scenarios meeting the validity criterion (n_overflows <= x),
	# appended to the auto-generated title (a plotly-safe place -- ggplotly
	# drops ggplot subtitles).
	valid_pct <- round(100 * mean(
	simulation_results_optimisation$n_overflows <= x_int, na.rm = TRUE))
	share_txt <- switch(lang,
	de = paste0(valid_pct, " % mit <= ", x_int, " Überläufen"),
	en = paste0(valid_pct, " % with <= ", x_int, " overflows"))
	if (is.null(title)) {
	title <- paste0(def_title, " (", share_txt, ")")
	if (y_var == "cost_per_evap_pct") {
	# name the reference of the marginal measure in the title
	title <- paste0(title, "\n", switch(lang,
	de = paste0("Referenz: minimale Evapotranspiration der gültigen ",
	"Szenarien (<= ", x_int, " Überläufe): ",
	round(evap_min, 1), " % (Szenario ",
	evap_min_scenario, ")"),
	en = paste0("Reference: minimum evapotranspiration of the valid ",
	"scenarios (<= ", x_int, " overflows): ",
	round(evap_min, 1), " % (scenario ",
	evap_min_scenario, ")")))
	}
	}

[claude]

Only cv$time[2] is updated — silently corrupts any Growth/Shading curve with more than 2 points.

The extension logic assumes the shipped Growth_1 / Shading_1 curves are exactly two-point ramps (start → end). For any curve with 3+ knots, points 3..n retain their original time values while point 2 is jumped to the rain end, producing a non-monotone (or otherwise physically wrong) curve that the engine reads without complaint.

Given the site templates the vignettes exercise are all 2-point ramps this hasn't bitten yet, but the length(cv$time) >= 2 guard reads as if any length were handled. Either:

• tighten the guard to length(cv$time) == 2 and stop/skip otherwise, so the assumption is enforced; or
• rescale all inner knots proportionally, e.g. cv$time <- cv$time * (max(timeseries_rain$time) / cv$time[length(cv$time)]).

At minimum, add a comment stating the 2-point assumption.

[claude]

warn.1d.NelderMead = FALSE does nothing here — the search is 3-D.

That control is only consulted when par has length 1 (see ?stats::optim), so on a 3-parameter starting point it is silently ignored. Not a bug, but the flag reads as if it were suppressing something and future-maintainers will wonder what. Suggest removing it, or replacing with a brief comment noting that the 1-D branch is intentionally impossible here.

Suggested change

	warn.1d.NelderMead = FALSE))
	stats::optim(starts[[si]], objective, method = "Nelder-Mead",
	control = list(maxit = 200, reltol = 1e-4))

[claude]

Wobble ladder only searches a geometric offset sequence (1, 2, 4, 8, …) below the upper edge — a strictly-monotone-but-decreasing-in-large-steps violation between the edge and lower bounds can be missed.

In discrete mode with axis_hi = length(grid) = 4 and wobble = 1, the ladder tests indices 3 then 2 (offset 1 then 2), then p = 0 <= axis_lo breaks. If only index 1 is feasible while 2 and 3 are not, the ladder declares infeasible even though a feasible level exists.

This is documented as a monotonicity assumption (the volume referee catches most violations) and the description at line 12–23 spells this out, so it's likely intentional. But the ladder's exponential search means the tolerance window shrinks the further below axis_hi the feasible anchor sits — worth calling out in the docstring so users know the wobble is only a neighbour fix, not a global rescue. Consider adding one sentence to the doc:

The ladder walks downward in offsets axis_tol, 2·axis_tol, 4·axis_tol, …; feasible anchors more than a few tolerance steps below the edge will only be found if strict monotonicity holds up to that step size.

[claude]

cost_per_evap_pct uses .data$cost_total in the numerator regardless of the active y_var.

The definition is cost_total / (evap - evap_min), but when the caller is going to filter by y_col = "cost_per_evap_pct" anyway, this is fine. Just worth verifying: is cost_total the intended numerator, or should it be the raw y_col value? If y_col == "cost_per_evap_pct" gets fed back through this expression on a subsequent invocation, the divide-by-evap would compound.

In current usage the numerator is cost_total unconditionally, so this is not a bug today — flagging as a defensive-review concern only. Consider a short comment: # always cost_total (not y_col) — the *derivation* of the marginal cost measure.

[claude]

u1 = (area - area_bounds[1]) / diff(area_bounds) — silently NaN if degenerate area_bounds = c(v, v) are ever passed.

Unlike u2 (guarded at line 423 with up > height_bounds[1]) and u3 (guarded for discrete via index; guarded for continuous via gb[2] > gb[1]), there's no diff(area_bounds) > 0 guard on the area axis. Passing area_bounds = c(100, 100) would produce u1 = NaN, which then propagates through pmin/pmax (both return NaN) and into stats::optim() where it errors with an opaque message.

Not encountered in the shipped defaults, and the stopifnot() at line 255 would arguably be the right place to add diff(area_bounds) > 0 if you want to fail loudly rather than by silent NaN.

[claude]

Nice reference-optimum design — worth extending one edge that isn't currently covered.

The synthetic-model tests are impressively thorough. One case that isn't exercised in either test-optimise_swale_design.R or test-optimise_swale_design_simultaneous.R: the volume referee firing end-to-end (a monotonicity violation flowing all the way to monotonicity_warning = TRUE in the output tibble). test-find_min_feasible.R covers it for the primitive, but the plumbing that assigns res_a$monotonicity_violation and res_s$monotonicity_violation and OR's them into mono_warn in solve_cell() (R/optimise_swale_design.R:268, 287, 304) is not tested — a synthetic run_fn that returns a real non-monotonicity (n_overflows AND overflow volume increase with area) would catch a regression in that wiring.

Non-blocking; just a suggestion for a follow-up commit.

[claude]

which(simulation_names == s_name) inside the lapply returns multiple indices if simulation_names has duplicates, then sprintf errors out.

Purely defensive: shipping code and tests never pass duplicate scenario names, so this only matters for hand-written calls. A simple fix that stays byte-identical for the common (unique) case:

messageText = sprintf("(%d/%d)) Reading results files for model run %s",
                      match(s_name, simulation_names),
                      length(simulation_names),
                      paths$dir_target_output),

match() returns just the first index — always length 1.

[claude]

size_cap semantics are only correct for size_by = "overflow_volume" — they invert the "larger = better" meaning for size_by = "evapotranspiration".

The clip-to-valid-region logic is justified for overflow volume: invalid scenarios have very large volumes that would compress the whole scale, so capping at the valid-region max keeps the valid points readable.

But for size_by = "evapotranspiration" (line 60–63 of the docstring: "larger points then mean more evapotranspiration, which is desirable"), invalid scenarios can have evapotranspiration higher or lower than valid ones — and if higher, capping them makes the desirable outliers look like the valid-region maximum, hiding the "outperformer that overflows too often" signal that the size scale is meant to convey.

Suggested guard: only apply the cap when size_by == "overflow_volume".

if (size_by == "overflow_volume") {
  valid_size <- df[[size_col]][!is.na(df$n_overflows) & df$n_overflows <= x_int]
  size_cap <- suppressWarnings(max(valid_size[is.finite(valid_size)]))
  if (!is.finite(size_cap) || size_cap <= 0) {
    size_cap <- suppressWarnings(max(df[[size_col]], na.rm = TRUE))
  }
  if (!is.finite(size_cap) || size_cap <= 0) size_cap <- 1
  df$size_plot <- pmin(df[[size_col]], size_cap)
} else {
  df$size_plot <- df[[size_col]]
}

[claude]

Sibling default mismatch: x = 5 here vs. x = 1 in plot_cost_vs_overflow_volume / plot_cost_vs_evaporation.

The docstring (line 22–27) advertises the tooltip as identical to plot_cost_vs_overflow_volume, and both plots derive evap_min from "scenarios that satisfy n_overflows <= x", so the same tooltip line will actually differ between the two plots at their defaults. If a user renders both in a vignette without overriding x, they'll see the same scenario reported with different cost_per_evap numbers side by side.

Either align the defaults (both to x = 1 or both to x = 5), or add a note to both docstrings that the derived tooltip line depends on x and users should pass the same value to both.

[claude]

Two more findings for this PR that I couldn't attach as inline comments because the target lines fall outside this PR's diff hunks (they're in unchanged prose / defaults of a heavily-changed file):

HIGH — R/plot_valid_design_space.R:164-172 — English title uses raw column names, no English axis-label dictionary

The English branch of txt$title (line 165-172) interpolates the raw x and y arguments (the column names), whereas the German branch above uses the resolved lab_x / lab_y. Combined with axis_labels_de (line 135-140) being German-only, plot_valid_design_space(..., lang = "en", x = "mulde_area", y = "mulde_height") produces:

"Valid solutions (Number of overflow events <= 5)\nin design space: mulde_area × mulde_height"

Suggested fix: route both languages through the shared default_param_labels(lang) dictionary and use lab_x / lab_y in the English paste0 too. Bonus: default_param_labels() covers filter_height, connected_area, rain_factor, lai, storage_type, bottom_hydraulicconductivity, all currently untranslated even in German.

HIGH — R/plot_wb_tradeoff_overflows.R:6-14, 38-42 — docstring off-by-one

The description says:

red for n_overflows >= x, shown as category ">=x"

and filter_n_gtx param doc:

If TRUE, scenarios with n_overflows >= x are removed before plotting.

But the code uses strict > (line 176 hi_lab <- paste0(">", x_int), line 182 n_overflows > x_int ~ hi_lab, line 175 n_overflows <= x_int), matching the sibling cost plots. A user setting x = 3 per the docs to drop >= 3 overflows actually keeps 0, 1, 2, 3. Since behaviour matches the siblings, the fix is docs-only: >=x → >x, ">=x" → ">x", and the filter_n_gtx param description.

---

These are the last two of a batch of 14 review findings on this PR; the other 12 posted successfully as inline comments.