Adapt ddwAgeCa to several iterations

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: FLBEIA
 Title: Bio-Economic Impact Assessment of Management Strategies using FLR
-Version: 1.16.1.16
-Date: 2024-06-28
+Version: 1.16.1.19
+Date: 2024-07-16
 Authors@R: c(person("FLBEIA", "Team", email = "flbeia@azti.es", role = c("aut","cre")))
 Description: A simulation toolbox that describes a fishery system under
   a Management Strategy Estrategy approach. The objective of the model is

NEWS

@@ -1,3 +1,18 @@
+#~~~~~~~~~~~~~~~~
+# 2024/07/12
+#~~~~~~~~~~~~~~~~
+FLBEIA 1.16.1.19
+- XXXXX
+
+
+
+#~~~~~~~~~~~~~~~~
+# 2024/07/09
+#~~~~~~~~~~~~~~~~
+FLBEIA 1.16.1.17
+- New price functions available (elasticPriceAge) and previous ones (elasticPrice) include alternative options now.
+
+
 #~~~~~~~~~~~~~~~~
 # 2024/02/20
 #~~~~~~~~~~~~~~~~

R/OM_1a2_DensityDependent_weight_functions.R

@@ -23,16 +23,21 @@ ddwAgeCa <- function(biol, stknm, year, season, ctrl, covars, ...) {
 
   require(nleqslv)
 
-  pars  <- ctrl$params[,,year,season,] # array[npar,nage,nyr,ns,nit]
+  na <- dim(biol@wt)[1]
+  it <- dim(biol@wt)[6]
+
+  pars  <- ctrl$params[,,year,season,] # array[npar,nage,1,1,nit]
   covnm <- ctrl$covnm                  # name of the covariate
 
-  wt.ref <- biol@wt[,year,,season,drop=TRUE]
-  ssna   <- (n(biol) * mat(biol) * exp(-spwn(biol) * m(biol)))[,year,,season,drop=TRUE] # exclude age0 (inmature)
-  sst    <- covars[[covnm]][stknm,year,,season,drop=TRUE]
+  wt.ref <- biol@wt[,year,,season,drop=TRUE] # [na,nit] or [na]
+  ssna   <- (n(biol) * mat(biol) * exp(-spwn(biol) * m(biol)))[,year,,season,drop=TRUE] # [na,nit] or [na] # to exclude age0 (immature)
+  sst    <- covars[[covnm]][stknm,year,,season,drop=TRUE]  # [nit] 
+
+  wage <- wt.ref * NA
 
   # ssb.ref <- quantSums(ssna * wt.ref[-1])
 
-  wfun <- function(x) {
+  wfun <- function(x, param, ssna) {
 
     # x[i] is log(wage_i)
     # x[i] = a_i + b_i * log(ssb) + c_i * sst
@@ -44,33 +49,59 @@ ddwAgeCa <- function(biol, stknm, year, season, ctrl, covars, ...) {
     logssb <- log(ssna[1] * exp(x[1]) + ssna[2] * exp(x[2]) +
                     ssna[3] * exp(x[3]) + ssna[4] * exp(x[4]))
 
-    y[1] <- x[1] - pars["b",1] * logssb - pars["a",1] - pars["c",1] * sst
-    y[2] <- x[2] - pars["b",2] * logssb - pars["a",2] - pars["c",2] * sst
-    y[3] <- x[3] - pars["b",3] * logssb - pars["a",3] - pars["c",3] * sst
-    y[4] <- x[4] - pars["b",4] * logssb - pars["a",4] - pars["c",4] * sst
+    y[1] <- x[1] - param["b",1] * logssb - param["a",1] - param["c",1] * sst
+    y[2] <- x[2] - param["b",2] * logssb - param["a",2] - param["c",2] * sst
+    y[3] <- x[3] - param["b",3] * logssb - param["a",3] - param["c",3] * sst
+    y[4] <- x[4] - param["b",4] * logssb - param["a",4] - param["c",4] * sst
 
     return(y)
   }
 
-  lwest <- nleqslv(log(wt.ref), fn = wfun, control = list(btol=1e-08, delta="newton"))
-
-  if (lwest$termcd != 1)
-    print(paste0(stknm," weights calculation message: ",lwest3$message))
-
-  if (any(wfun(lwest$x)>1e-08))
-    stop(paste("Issues with mean weight calculation for", stknm))
-
-  wage <- exp(lwest$x)
-
-  # check
-  logssb <- sum(ssna * wage)
-  di <- numeric(4)
-  for (i in 1:4) {
-    di[i] <- wage[i] - exp(pars["a",i] + pars["b",i] * logssb + pars["c",i] * sst)
+  if (it == 1) {
+
+    lwest <- nleqslv(log(wt.ref), fn = wfun, param=pars, ssna=ssna, control = list(btol=1e-08, delta="newton"))
+
+    if (lwest$termcd != 1)
+      print(paste0(stknm," weights calculation message: ",lwest3$message))
+
+    if (any(wfun(lwest$x, param=pars, ssna=ssna)>1e-08))
+      stop(paste("Issues with mean weight calculation for", stknm))
+
+    wage <- exp(lwest$x)
+
+    # # check
+    # logssb <- sum(ssna * wage)
+    # di <- numeric(na)
+    # for (a in 1:na)
+    #   di[a] <- wage[a] - exp(pars["a",a] + pars["b",a] * logssb + pars["c",a] * sst)
+    # if (any(round(di-wage,4) != 0))
+    #   stop(paste("Issues with density-dependent weights-at-age for",stknm))
+
+  } else {
+
+    for (i in 1:it) {
+
+      lwest <- nleqslv(log(wt.ref[,i]), fn = wfun, param=pars[,,i], ssna=ssna[,i], control = list(btol=1e-08, delta="newton"))
+
+      if (lwest$termcd != 1)
+        print(paste0(stknm," weights calculation message: ",lwest3$message))
+
+      if (any(wfun(lwest$x, param=pars[,,i], ssna=ssna[,i])>1e-08))
+        stop(paste("Issues with mean weight calculation for", stknm))
+
+      wage[,i] <- exp(lwest$x)
+    }
+
+    # # check
+    # logssb <- apply(ssna * wage, 2, sum)
+    # di <- array(NA, dim=c(na,it))
+    # for (a in 1:na)
+    #   di[a,] <- wage[a,] - exp(pars["a",a,] + pars["b",a,] * logssb + pars["c",a,] * sst)
+    # if (any(round(di-wage,4) != 0))
+    #   stop(paste("Issues with density-dependent weights-at-age for",stknm))
+
   }
-  if (any(round(di-wage,8) != 0))
-    stop(paste("Issues with density-dependent weights-at-age for",stknm))
-
+
   # wt change
   wt.chg <- wage/wt.ref
 

README.md

@@ -1,6 +1,6 @@
 # FLBEIA
-- Version: 1.16.1.18
-- Date: 2024-07-09
+- Version: 1.16.1.19
+- Date: 2024-07-16
 - Author: Dorleta GARCIA <dgarcia@azti.es>; FLBEIA Team <flbeia@azti.es>
 - Maintainer: Dorleta GARCIA, AZTI & FLBEIA Team
 - Repository: <https://github.com/flr/FLBEIA/>