R Under development (unstable) (2024-03-18 r86148 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> if (!interactive()) options(warn=2, error = function() { sink(stderr()) ; traceback(3) ; q(status = 1) })
> library(magrittr)
> library(unittest)
> 
> library(gadget3)
> 
> ok_group('g3a_mature_constant', {
+     cmp_code <- function(a, b) ut_cmp_identical(rlang::f_rhs(a), rlang::f_rhs(b))
+ 
+     ok(cmp_code(
+         g3a_mature_constant(),
+         ~1/(1 + exp(0)) ), "No arguments")
+ 
+     # Invalid combinations
+     ok(ut_cmp_error(g3a_mature_constant(alpha = 2), "l50"), "Missing l50")
+     ok(ut_cmp_error(g3a_mature_constant(beta = 2), "a50"), "Missing a50")
+     ok(ut_cmp_error(g3a_mature_constant(gamma = 2), "k50"), "Missing k50")
+ 
+     # Single alpha/beta/gamma
+     ok(cmp_code(
+         g3a_mature_constant(alpha = 28, l50 = 24),
+         ~1/(1 + exp((0 - (28) * (stock__midlen - (24))))) ), "alpha = 28, l50 = 24")
+     ok(cmp_code(
+         g3a_mature_constant(beta = 18, a50 = 83),
+         ~1/(1 + exp((0 - (18) * (age - (83))))) ), "beta = 18, a50 = 83")
+     ok(cmp_code(
+         g3a_mature_constant(gamma = 82, k50 = 27),
+         ~1/(1 + exp((0 - (82) * (stock_ss(stock__wgt) - (27))))) ), "gamma = 82, k50 = 27")
+ 
+     # Zero will also remove conditions
+     ok(cmp_code(
+         g3a_mature_constant(alpha = 28, l50 = 24, beta = 0, gamma = 0),
+         ~1/(1 + exp((0 - (28) * (stock__midlen - (24))))) ), "alpha = 28, l50 = 24")
+     ok(cmp_code(
+         g3a_mature_constant(beta = 18, a50 = 83, alpha = 0, gamma = 0),
+         ~1/(1 + exp((0 - (18) * (age - (83))))) ), "beta = 18, a50 = 83")
+ 
+     # Can combine parameters
+     ok(cmp_code(
+         g3a_mature_constant(alpha = 73, l50 = 16, beta = 54, a50 = 32),
+         ~1/(1 + exp(((0 - (73) * (stock__midlen - (16)))
+                       - (54) * (age - (32)))
+                       ))), "alpha = 73, l50 = 16, beta = 54, a50 = 32")
+     ok(cmp_code(
+         g3a_mature_constant(beta = 73, a50 = 67, gamma = 39, k50 = 73),
+         ~1/(1 + exp(((0 - (73) * (age - (67)))
+                       - (39) * (stock_ss(stock__wgt) - (73))
+                       )))), "beta = 73, a50 = 67, gamma = 39, k50 = 73")
+ 
+     # alpha/beta/gamma can be formula
+     ok(cmp_code(
+         g3a_mature_constant(alpha = ~g3_param("ling.mat1"), l50 = ~g3_param("ling.mat2")),
+         ~1/(1 + exp((0 - (g3_param("ling.mat1")) * (stock__midlen - (g3_param("ling.mat2"))))))), "alpha = formula")
+ 
+     ok(cmp_code(
+         g3a_mature_constant(beta = ~g3_param("ling.mat1"), a50 = ~g3_param("ling.mat2")),
+         ~1/(1 + exp((0 - (g3_param("ling.mat1")) * (age - (g3_param("ling.mat2"))))))), "beta = formula")
+ 
+     ok(cmp_code(
+         g3a_mature_constant(gamma = ~g3_param("ling.mat1"), k50 = ~g3_param("ling.mat2")),
+         ~1/(1 + exp((0 - (g3_param("ling.mat1")) * (stock_ss(stock__wgt) - (g3_param("ling.mat2"))))))), "gamma = formula")
+ })
# g3a_mature_constant
ok - No arguments
ok - Missing l50
ok - Missing a50
ok - Missing k50
ok - alpha = 28, l50 = 24
ok - beta = 18, a50 = 83
ok - gamma = 82, k50 = 27
ok - alpha = 28, l50 = 24
ok - beta = 18, a50 = 83
ok - alpha = 73, l50 = 16, beta = 54, a50 = 32
ok - beta = 73, a50 = 67, gamma = 39, k50 = 73
ok - alpha = formula
ok - beta = formula
ok - gamma = formula
> 
> ok_group('g3a_mature', {
+     stock_imm <- g3_stock('stock_imm', seq(10, 40, 10)) %>% g3s_age(3, 7)
+     stock_mat1 <- g3_stock('stock_mat1', seq(10, 40, 10)) %>% g3s_age(4, 7)
+     stock_mat2 <- g3_stock('stock_mat2', seq(10, 40, 10)) %>% g3s_age(5, 7)
+ 
+     ok(ut_cmp_error(
+         g3a_mature(stock_imm, ~g3_param_vector("maturity"), list(stock_mat1, stock_mat2), output_ratios = c(9,9,9)),
+         "output_ratios"), "Length of output_ratios must match")
+     ok(ut_cmp_error(
+         g3a_mature(stock_imm, ~g3_param_vector("maturity"), list(stock_mat1, stock_mat2), output_ratios = c(9,9)),
+         "output_ratios"), "output_ratios must sum to 1")
+ 
+     actions <- list(
+         g3a_time(2000, 2000, project_years = 0),
+         g3a_initialconditions(stock_imm, ~max(6L - age, 1L) * g3_param_vector("imm_init_num"), ~g3_param_vector("imm_init_wgt")),
+         g3a_initialconditions(stock_mat1, ~max(6L - age, 1L) * g3_param_vector("mat1_init_num"), ~g3_param_vector("mat1_init_wgt")),
+         g3a_initialconditions(stock_mat2, ~max(6L - age, 1L) * g3_param_vector("mat2_init_num"), ~g3_param_vector("mat2_init_wgt")),
+         g3a_mature(stock_imm,
+             ~g3_param_vector("maturity"),
+             list(stock_mat1, stock_mat2),
+             output_ratios = list(~g3_param("ratio_mat1"), ~g3_param("ratio_mat2")),
+             run_f = ~g3_param("run_f") == 1),
+         list(
+             '999' = ~{
+                 REPORT(stock_imm__num)
+                 REPORT(stock_imm__wgt)
+                 REPORT(stock_mat1__num)
+                 REPORT(stock_mat1__wgt)
+                 REPORT(stock_mat2__num)
+                 REPORT(stock_mat2__wgt)
+                 nll <- nll + g3_param('x')
+                 return(nll)
+             }))
+ 
+     # Compile model
+     model_fn <- g3_to_r(actions, trace = FALSE)
+     # model_fn <- edit(model_fn)
+     params <- attr(model_fn, 'parameter_template')
+     params$imm_init_num <- c(101, 102, 103, 104)
+     params$imm_init_wgt <- c(1, 2, 3, 4)
+     params$maturity <- c(0, 0, 1, 0)
+     params$mat1_init_num <- c(10, 20, 30, 40)
+     params$mat1_init_wgt <- c(9, 8, 7, 6)
+     params$mat2_init_num <- c(0, 0, 0, 0)
+     params$mat2_init_wgt <- c(0, 0, 0, 0)
+     params$ratio_mat1 <- 0.5
+     params$ratio_mat2 <- 0.5
+     params$run_f <- 1  # TODO: if(param) is naughty
+ 
+     if (nzchar(Sys.getenv('G3_TEST_TMB'))) {
+         model_cpp <- g3_to_tmb(actions, trace = FALSE)
+         # model_cpp <- edit(model_cpp)
+         model_tmb <- g3_tmb_adfun(model_cpp, params, compile_flags = c("-O0", "-g"))
+     } else {
+         writeLines("# skip: not compiling TMB model")
+     }
+ 
+     ok_group("Move all of length 30 into mat1/mat2", {
+         params <- attr(model_fn, 'parameter_template')
+         params$imm_init_num <- c(101, 102, 103, 104)
+         params$imm_init_wgt <- c(1, 2, 3, 4)
+         params$maturity <- c(0, 0, 1, 0)
+         params$mat1_init_num <- c(10, 20, 30, 40)
+         params$mat1_init_wgt <- c(9, 8, 7, 6)
+         params$mat2_init_num <- c(0, 0, 0, 0)
+         params$mat2_init_wgt <- c(0, 0, 0, 0)
+         params$ratio_mat1 <- 0.5
+         params$ratio_mat2 <- 0.5
+         params$run_f <- 1
+ 
+         result <- model_fn(params)
+         r <- attributes(result)
+         ok(ut_cmp_equal(
+             as.vector(r$stock_imm__num[,'age5']),
+             c(101, 102, 0, 104)), "stock_imm__num")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_imm__num[length = '30:40',]),
+             c(103 * 3, (103 * 2) / 2, 0, 0, 0)), "stock_imm__num: age3 and half of age4 left behind, since they don't fit in mature stocks")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat1__num[,'age5']),
+             c(10, 20, 30 + (103 / 2), 40)), "stock_mat1__num (NB: Added to existing numbers)")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat1__num[length = '30:40',]),
+             c(60 + (103 * 2) / 2, 81.5, 81.5, 81.5)), "stock_mat1__num: Got half of age4")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat2__num[,'age5']),
+             c(0, 0, 103 / 2, 0)), "stock_mat2__num")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat2__num[length = '30:40',]),
+             c(51.5, 51.5, 51.5)), "stock_mat2__num: Just age 5/6/7")
+ 
+         if (nzchar(Sys.getenv('G3_TEST_TMB'))) {
+             param_template <- attr(model_cpp, "parameter_template")
+             param_template$value <- params[param_template$switch]
+             gadget3:::ut_tmb_r_compare(model_fn, model_tmb, param_template)
+         }
+     })
+ 
+     ok_group("Move all of length 30 into 90% mat1, 10% mat2", {
+         params <- attr(model_fn, 'parameter_template')
+         params$imm_init_num <- c(101, 102, 103, 104)
+         params$imm_init_wgt <- c(1, 2, 3, 4)
+         params$maturity <- c(0, 0, 1, 0)
+         params$mat1_init_num <- c(10, 20, 30, 40)
+         params$mat1_init_wgt <- c(9, 8, 7, 6)
+         params$mat2_init_num <- c(0, 0, 0, 0)
+         params$mat2_init_wgt <- c(0, 0, 0, 0)
+         params$ratio_mat1 <- 0.9
+         params$ratio_mat2 <- 0.1
+         params$run_f <- 1
+ 
+         result <- model_fn(params)
+         r <- attributes(result)
+         ok(ut_cmp_equal(
+             as.vector(r$stock_imm__num[,'age5']),
+             c(101, 102, 0, 104)), "stock_imm__num")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat1__num[,'age5']),
+             c(10, 20, 30 + (103 * 0.9), 40)), "stock_mat1__num (NB: Added to existing numbers)")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat2__num[,'age5']),
+             c(0, 0, 103 * 0.1, 0)), "stock_mat2__num")
+ 
+         if (nzchar(Sys.getenv('G3_TEST_TMB'))) {
+             param_template <- attr(model_cpp, "parameter_template")
+             param_template$value <- params[param_template$switch]
+             gadget3:::ut_tmb_r_compare(model_fn, model_tmb, param_template)
+         }
+     })
+ 
+     ok_group("Move 50% of length 30, 75% of length 40", {
+         params <- attr(model_fn, 'parameter_template')
+         params$imm_init_num <- c(101, 102, 103, 104)
+         params$imm_init_wgt <- c(1, 2, 3, 4)
+         params$maturity <- c(0, 0, 0.5, 0.75)
+         params$mat1_init_num <- c(10, 20, 30, 40)
+         params$mat1_init_wgt <- c(9, 8, 7, 6)
+         params$mat2_init_num <- c(0, 0, 0, 0)
+         params$mat2_init_wgt <- c(0, 0, 0, 0)
+         params$ratio_mat1 <- 0.5
+         params$ratio_mat2 <- 0.5
+         params$run_f <- 1
+ 
+         result <- model_fn(params)
+         r <- attributes(result)
+         ok(ut_cmp_equal(
+             as.vector(r$stock_imm__num[,'age5']),
+             c(101, 102, 103 * 0.5, 104 * 0.25)), "stock_imm__num")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat1__num[,'age5']),
+             c(10, 20, 30 + (103 * 0.5 / 2), 40 + (104 * 0.75 / 2))), "stock_mat1__num (NB: Added to existing numbers)")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat2__num[,'age5']),
+             c(0, 0, 103 * 0.5 / 2, 104 * 0.75 / 2)), "stock_mat2__num")
+ 
+         if (nzchar(Sys.getenv('G3_TEST_TMB'))) {
+             param_template <- attr(model_cpp, "parameter_template")
+             param_template$value <- params[param_template$switch]
+             gadget3:::ut_tmb_r_compare(model_fn, model_tmb, param_template)
+         }
+     })
+ 
+     ok_group("Disable with run_f = 0", {
+         params <- attr(model_fn, 'parameter_template')
+         params$imm_init_num <- c(101, 102, 103, 104)
+         params$imm_init_wgt <- c(1, 2, 3, 4)
+         params$maturity <- c(0, 0, 0.5, 0.75)
+         params$mat1_init_num <- c(10, 20, 30, 40)
+         params$mat1_init_wgt <- c(9, 8, 7, 6)
+         params$mat2_init_num <- c(0, 0, 0, 0)
+         params$mat2_init_wgt <- c(0, 0, 0, 0)
+         params$ratio_mat1 <- 0.5
+         params$ratio_mat2 <- 0.5
+         params$run_f <- 0  # NB: Turned off
+ 
+         result <- model_fn(params)
+         r <- attributes(result)
+         ok(ut_cmp_equal(
+             as.vector(r$stock_imm__num[,'age5']),
+             c(101, 102, 103, 104)), "stock_imm__num same as start")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat1__num[,'age5']),
+             c(10, 20, 30, 40)), "stock_mat1__num same as start")
+         ok(ut_cmp_equal(
+             as.vector(r$stock_mat2__num[,'age5']),
+             c(0, 0, 0, 0)), "stock_mat2__num same as start")
+ 
+         if (nzchar(Sys.getenv('G3_TEST_TMB'))) {
+             param_template <- attr(model_cpp, "parameter_template")
+             param_template$value <- params[param_template$switch]
+             gadget3:::ut_tmb_r_compare(model_fn, model_tmb, param_template)
+         }
+     })
+ })
# g3a_mature
ok - Length of output_ratios must match
ok - output_ratios must sum to 1
# skip: not compiling TMB model
# Move all of length 30 into mat1/mat2
ok - stock_imm__num
ok - stock_imm__num: age3 and half of age4 left behind, since they don't fit in mature stocks
ok - stock_mat1__num (NB: Added to existing numbers)
ok - stock_mat1__num: Got half of age4
ok - stock_mat2__num
ok - stock_mat2__num: Just age 5/6/7
# Move all of length 30 into 90% mat1, 10% mat2
ok - stock_imm__num
ok - stock_mat1__num (NB: Added to existing numbers)
ok - stock_mat2__num
# Move 50% of length 30, 75% of length 40
ok - stock_imm__num
ok - stock_mat1__num (NB: Added to existing numbers)
ok - stock_mat2__num
# Disable with run_f = 0
ok - stock_imm__num same as start
ok - stock_mat1__num same as start
ok - stock_mat2__num same as start
> 
> ok_group('movement with non-conforming stocks', {
+     # Transform required from a -> b+c, make sure we can form a model
+     prey_a <- g3_stock('prey_a', seq(20, 160, 4)) %>% g3s_age(11, 12)
+     prey_b <- g3_stock('prey_b', seq(20, 160, 4)) %>% g3s_age(11, 15)
+     prey_c <- g3_stock('prey_c', seq(20, 164, 4)) %>% g3s_age(11, 15)
+     # Build model, if missing vars are present the warning will be translated into error
+     m <- g3_to_r(list(
+         g3a_time(2000, 2002, step_lengths = c(6, 6), project_years = 0),
+         g3a_age(prey_a, output_stocks = list(prey_b, prey_c), output_ratios = c(0.5,0.5))))
+ })
# movement with non-conforming stocks
> 
> proc.time()
   user  system elapsed 
   1.03    0.09    1.10 
# Looks like you passed all 31 tests.