Add module3: cost of storing fertilizers

DESCRIPTION

@@ -33,6 +33,7 @@ Imports:
     R6,
     torch
 Suggests: 
+    devtools,
     testthat (>= 3.0.0)
 Config/testthat/edition: 3
 Depends: 

NEWS.md

@@ -2,5 +2,5 @@
 
 ## Added
 * Adds `apus` object with the functions `addField`, `trainModel` and `optimizeFertilizerChoice`
-* Adds for the cost function module 1 and 2
+* Adds for the cost function module 1,5 and 6
 * Adds default table for `cultivations`, `fertilizers`, `parameters` and `fines`

R/apus.R

@@ -133,6 +133,48 @@ Apus <- R6::R6Class(
       return(TRUE)
     },
 
+    #' @description
+    #' Update a fertilizer
+    #'
+    #' @param p_id (character) ID of the fertilizer
+    #' @param p_price (number)
+    #' @param p_stored (number)
+    #' @param p_storage_available (number)
+    #'
+    #' @export
+    updateFertilizer = function(p_id, p_price = NA_real_, p_stored = NA_real_, p_storage_available = NA_real_) {
+
+      # Check arguments ---------------------------------------------------------
+      #TODO
+
+      # update the fertilizers table with data for that fertilizer --------------
+      fertilizers.new <- self$fertilizers
+      p_id.new <- p_id
+
+      # Update price
+      if(! is.na(p_price)) {
+        p_price.new <- p_price
+        fertilizers.new[p_id == p_id.new, p_price := p_price.new]
+      }
+
+      # Update stored amount of fertilizer
+      if(! is.na(p_stored)) {
+        p_stored.new <- p_stored
+        fertilizers.new[p_id == p_id.new, p_stored := p_stored.new]
+      }
+
+      # Update available storages for fertilizer
+      if(! is.na(p_storage_available)) {
+        p_storage_available.new <- p_storage_available
+        fertilizers.new[p_id == p_id.new, p_storage_available := p_storage_available.new]
+      }
+
+      # Return back updated fertilizer
+      self$fertilizers <- fertilizers.new
+
+      return(TRUE)
+    },
+
     #' @description
     #' Train a model to
     #'
@@ -168,6 +210,7 @@ Apus <- R6::R6Class(
       dataset.train <- createApusDataset(farms = NULL,
                                          cultivation = self$cultivation,
                                          fertilizers = self$fertilizers,
+                                         fines = self$fines,
                                          fields_max = self$fields_max,
                                          device = device)
 
@@ -177,6 +220,7 @@ Apus <- R6::R6Class(
       dataset.valid <- createApusDataset(farms = farms.valid,
                                          cultivation = self$cultivation,
                                          fertilizers = self$fertilizers,
+                                         fines = self$fines,
                                          fields_max = self$fields_max,
                                          device = device)
 
@@ -217,7 +261,7 @@ Apus <- R6::R6Class(
         fields[is.na(fields)] <- 0
       }
 
-      dataset <- createApusDataset(farms = fields, cultivations = self$cultivations, fertilizers = self$fertilizers, fields_max = self$fields_max, device = self$device)
+      dataset <- createApusDataset(farms = fields, cultivations = self$cultivations, fertilizers = self$fertilizers, fines = self$fines, fields_max = self$fields_max, device = self$device)
       dl <- torch::dataloader(dataset, batch_size = 1)
 
 

R/dataset.R

@@ -6,6 +6,7 @@
 #' @param farms (data.table)
 #' @param cultivations (data.table)
 #' @param fertilizers (data.table)
+#' @param fines (data.table)
 #' @param fields_max (integer)
 #' @param device (character)
 #'
@@ -15,7 +16,7 @@
 #' @import torch
 #'
 #'@export
-createApusDataset <- function(farms = NULL, cultivations, fertilizers, fields_max, device) {
+createApusDataset <- function(farms = NULL, cultivations, fertilizers, fines, fields_max, device) {
 
   transformfieldsToTensor = createSyntheticfields = code = fields_count = self = NULL
   size = value_max = value_min = p_price = p_stored = b_id_farm = b_id_field = NULL
@@ -27,7 +28,7 @@ createApusDataset <- function(farms = NULL, cultivations, fertilizers, fields_ma
   apus_dataset <- torch::dataset(
     name = "apus_dataset",
 
-    initialize = function(farms = NULL, cultivations, fertilizers, fields_max, device) {
+    initialize = function(farms = NULL, cultivations, fertilizers, fines, fields_max, device) {
 
       # Check arguments -----------------------------------------------------
       # TODO
@@ -45,12 +46,10 @@ createApusDataset <- function(farms = NULL, cultivations, fertilizers, fields_ma
         self$farms_count <- 100
       }
 
-      # Set temporary
-      fertilizers[, p_stored := 0]
-      fertilizers[, p_price := 1]
-
-      fertilizers <- fertilizers[, c('p_stored', 'p_price', 'p_n_rt', 'p_n_wc', 'p_p_rt', 'p_k_rt')]
+      fines <- dcast(fines, . ~ norm, value.var = 'fine')[, 2:4]
+      self$fines <- torch::torch_tensor(as.matrix(fines), device = device)
 
+      fertilizers <- fertilizers[, c('p_stored', 'p_price', 'p_n_rt', 'p_n_wc', 'p_p_rt', 'p_k_rt', 'p_type_manure', 'p_p_wcl', 'p_storage_cost', 'p_storage_capacity', 'p_storage_available')]
       self$fertilizers <- torch::torch_tensor(as.matrix(fertilizers), device = device)
     },
 
@@ -63,7 +62,7 @@ createApusDataset <- function(farms = NULL, cultivations, fertilizers, fields_ma
       }
       t.fields <- transformFieldsToTensor(farms, self$device)
 
-      return(list(fields = t.fields, fertilizers = self$fertilizers))
+      return(list(fields = t.fields, fertilizers = self$fertilizers, fines = self$fines))
     },
 
     .length = function() {
@@ -73,7 +72,7 @@ createApusDataset <- function(farms = NULL, cultivations, fertilizers, fields_ma
 
 
   # Create torch dataset for apus -------------------------------------------
-  dataset <- apus_dataset(farms = farms, cultivations = cultivations, fertilizers = fertilizers, fields_max = fields_max, device = device)
+  dataset <- apus_dataset(farms = farms, cultivations = cultivations, fertilizers = fertilizers, fines = fines, fields_max = fields_max, device = device)
 
   return(dataset)
 }

R/model.R

@@ -106,14 +106,14 @@ createApusModel <- function(dataset.train, dataset.valid, width = 12, layers = 1
 
       cli::cli_progress_bar(paste0('Training model [', epoch, '/', epochs, ']'), total = dl.train$.length())
 
-      # For testing
+      # For developing
       # b <- dl.train$.iter()
       # b <- b$.next()
 
       # Forward pass
       optimizer$zero_grad()
       doses <- model(b$fields, b$fertilizers)
-      cost <- calculateCost(doses, b$fields, b$fertilizers)
+      cost <- calculateCost(doses, b$fields, b$fertilizers, b$fines)
 
       # Backward pass
       cost$backward()
@@ -133,13 +133,13 @@ createApusModel <- function(dataset.train, dataset.valid, width = 12, layers = 1
 
       cli::cli_progress_bar(paste0('Validating model [', epoch, '/', epochs, ']'), total = dl.valid$.length())
 
-      # For testing
+      # For developing
       # b <- dl.valid$.iter()
       # b <- b$.next()
 
       # Forward pass
       doses <- model(b$fields, b$fertilizers)
-      cost <- calculateCost(doses, b$fields, b$fertilizers)
+      cost <- calculateCost(doses, b$fields, b$fertilizers, b$fines)
 
       losses.validation <- c(losses.validation, cost$item())
 
@@ -157,7 +157,7 @@ createApusModel <- function(dataset.train, dataset.valid, width = 12, layers = 1
   return(model)
 }
 
-calculateCost <- function(doses, fields, fertilizers, reduce_batches = TRUE) {
+calculateCost <- function(doses, fields, fertilizers, fines, reduce_batches = TRUE) {
 
 
   # Check arguments ---------------------------------------------------------
@@ -168,12 +168,20 @@ calculateCost <- function(doses, fields, fertilizers, reduce_batches = TRUE) {
   module1 <- calculateCostModule1(doses, fields, fertilizers)
 
 
-  # Module 4: Revenue from harvested crops ----------------------------------
-  module4 <- calculateRevenueModule4(doses, fields, fertilizers)
+  # Module 3: Cost of storing fertilizers -----------------------------------
+  module3 <- calculateCostModule3(doses, fields, fertilizers)
+
+
+  # Module 5: Revenue from harvested crops ----------------------------------
+  module5 <- calculateRevenueModule5(doses, fields, fertilizers)
+
+
+  # Module 6: Penalty for exceeding legal limit -----------------------------
+  module6 <- calculatePenaltyModule6(doses, fields, fertilizers, fines)
 
 
   # Combine the modules -----------------------------------------------------
-  cost <- torch::torch_zeros(dim(module1)) + module1 - module4
+  cost <- torch::torch_zeros(dim(module1)) + module1 + module3 - module5 + module6
 
 
   # Convert to € / ha -------------------------------------------------------
@@ -211,8 +219,31 @@ calculateCostModule1 <- function(doses, fields, fertilizers) {
   return(module1)
 }
 
-# Module 4: Revenue from harvested crops  ------------------------------------
-calculateRevenueModule4 <- function(doses, fields, fertilizers) {
+# Module 3: Cost of storing fertilizers -------------------------------------
+calculateCostModule3 <- function(doses, fields, fertilizers) {
+
+  # Sum dose per fertilizer -------------------------------------------------
+  fields.b_area <- torch::torch_unsqueeze(fields[,,1], -1)
+  fields.dose <- fields.b_area * doses
+  fertilizers.dose <- torch::torch_sum(doses, dim = 2L)
+
+
+  # Calculate requires storage places for fertilizer ------------------------
+  fertilizers.storage_capacity <- fertilizers[,,10]
+  fertilizers.storage_cost <- fertilizers[,,9]
+  fertilizers.storage_available <- fertilizers[,,11]
+  fertilizers.storages <- torch::torch_ceil(fertilizers.dose / fertilizers.storage_capacity)
+  fertilizers.cost <- (fertilizers.storages - fertilizers.storage_available) * fertilizers.storage_cost
+
+
+  # Sum cost for farm -------------------------------------------------------
+  module3 <- torch::torch_sum(fertilizers.cost, dim = 2L)
+
+  return(module3)
+}
+
+# Module 5: Revenue from harvested crops  ------------------------------------
+calculateRevenueModule5 <- function(doses, fields, fertilizers) {
 
   # Calculate N dose per fields
   fertilizers.p_n_rt <- fertilizers[,,3]
@@ -268,9 +299,71 @@ calculateRevenueModule4 <- function(doses, fields, fertilizers) {
   fields.b_area <- fields[,,1]
   fields.b_lu_yield <- fields[,,8]
   fields.b_lu_price <- fields[,,9]
-  module4 <- fields.b_area *  fields.b_lu_yield * fields.b_lu_price * fields.d_realized
-  module4 <- torch::torch_sum(module4, dim = 2L)
+  module5 <- fields.b_area *  fields.b_lu_yield * fields.b_lu_price * fields.d_realized
+  module5 <- torch::torch_sum(module5, dim = 2L)
+
+  return(module5)
+}
 
-  return(module4)
+# Module 6: Penalties in case of exceeding legal limits -----------------------
+calculatePenaltyModule6 <- function(doses, fields, fertilizers, fines) {
+
+  # Calculate d_n_norm_man per field
+  fertilizers.p_n_rt <- fertilizers[,,3]
+  fertilizers.p_type_manure <- fertilizers[,,7]
+  fertilizers.p_n_manure <- fertilizers.p_n_rt * fertilizers.p_type_manure
+  fertilizers.p_n_manure <- torch::torch_unsqueeze(fertilizers.p_n_manure, 2)
+  fertilizers.p_n_manure <- torch::torch_repeat_interleave(fertilizers.p_n_manure, repeats = dim(doses)[2], dim =2)
+  fields.fertilizers.dose.n_manure <- doses * fertilizers.p_n_manure
+  fields.dose.n_manure <- torch::torch_sum(fields.fertilizers.dose.n_manure, dim = 3)
+  farms.dose.n_manure <- torch::torch_sum(fields.dose.n_manure, dim = 2)
+
+  fields.d_n_norm_man <- fields[,,6]
+  fields.b_area <- fields[,,1]
+  fine.d_n_norm_man <- fines[,1,2]
+  farms.d_n_norm_man <- torch::torch_sum(fields.b_area * fields.d_n_norm_man, dim = 2L)
+  farms.exceeding.d_n_norm_man <- torch::torch_relu(farms.dose.n_manure - farms.d_n_norm_man)
+  farms.penalty.d_n_norm_man <- farms.exceeding.d_n_norm_man * fine.d_n_norm_man
+
+  # Calculate d_n_norm per field
+  fertilizers.p_n_rt <- fertilizers[,,3]
+  fertilizers.p_n_wc <- fertilizers[,,4] # TODO Replace with p_n_wcl
+  fertilizers.p_n_workable <- fertilizers.p_n_rt * fertilizers.p_n_wc
+  fertilizers.p_n_workable <- torch::torch_unsqueeze(fertilizers.p_n_workable, 2)
+  fertilizers.p_n_workable <- torch::torch_repeat_interleave(fertilizers.p_n_workable, repeats = dim(doses)[2], dim =2)
+  fields.fertilizers.dose.n_workable <- doses * fertilizers.p_n_workable
+  fields.dose.n_workable <- torch::torch_sum(fields.fertilizers.dose.n_workable, dim = 3)
+  farms.dose.n_workable <- torch::torch_sum(  fields.dose.n_workable , dim = 2)
+
+  fields.d_n_norm <- fields[,,5]
+  fields.b_area <- fields[,,1]
+  fine.d_n_norm <- fines[,1,1]
+  farms.d_n_norm <- torch::torch_sum(fields.b_area * fields.d_n_norm, dim = 2L)
+  farms.exceeding.d_n_norm <- torch::torch_relu(farms.dose.n_workable - farms.d_n_norm)
+  farms.penalty.d_n_norm <- farms.exceeding.d_n_norm * fine.d_n_norm
+
+  # Calculate d_p_norm per field
+  fertilizers.p_p_rt <- fertilizers[,,5]
+  fertilizers.p_p_wcl <- fertilizers[,,8]
+  fertilizers.p_p_legal <- fertilizers.p_p_rt * fertilizers.p_p_wcl
+  fertilizers.p_p_legal <- torch::torch_unsqueeze(fertilizers.p_p_legal, 2)
+  fertilizers.p_p_legal <- torch::torch_repeat_interleave(fertilizers.p_p_legal, repeats = dim(doses)[2], dim =2)
+  fields.fertilizers.dose.p_legal <- doses * fertilizers.p_p_legal
+  fields.dose.p_legal <- torch::torch_sum(fields.fertilizers.dose.p_legal, dim = 3)
+  farms.dose.p_legal <- torch::torch_sum(  fields.dose.p_legal, dim = 2)
+
+  fields.d_p_norm <- fields[,,7]
+  fields.b_area <- fields[,,1]
+  fine.d_p_norm <- fines[,1,3]
+  farms.d_p_norm <- torch::torch_sum(fields.b_area * fields.d_p_norm, dim = 2L)
+  farms.exceeding.d_p_norm <- torch::torch_relu(farms.dose.p_legal - farms.d_p_norm)
+  farms.penalty.d_p_norm <- farms.exceeding.d_p_norm * fine.d_p_norm
+
+
+  # Combine the penalties
+  module6 <- farms.penalty.d_n_norm_man + farms.penalty.d_n_norm + farms.penalty.d_p_norm
+
+  return(module6)
 }
 
+

README.md

@@ -57,6 +57,17 @@ apus$addField(
 )
 ```
 
+To update a fertilizer with specific values for your situation use the `updateFertilizer` function
+
+```
+apus$updateFertilizer(
+  p_id = 17,
+  p_price = -5,
+  p_stored = 15000,
+  p_storage_available = 1
+)
+```
+
 To get a fertilizer advice you need a model first. A model can be trained with with the function `trainModel()`
 
 `apus$trainModel()` 
@@ -73,13 +84,18 @@ For the v1 version of `apus` we plan to develop to following features:
 
 * [ ] Import and export trained models
 * [ ] Include a trained base model 
+* [x] Add function to train model
 * [ ] Enable fine-tuning of (base) models
-* [ ] Include cost function for module 2: Cost of storing fertilizers
-* [ ] Include cost function for module 3: Cost of applying fertilizers
-* [ ] Include cost function for module 5: Penalties in case of exceeding legal limits
-* [ ] Include cost function for module 6: Cost of greenhouse gas emissions
-* [ ] Include realistic cultivation response curves from module 4
-* [ ] Add other nutrients and organic matter to module 4
+* [x] Include cost function for module 1: Purchase of fertilizers
+* [ ] Include cost function for module 2: Disposal of manure
+* [ ] Include cost function for module 3: Cost of storing fertilizers
+* [ ] Include cost function for module 4: Cost of applying fertilizers
+* [x] Include cost function for module 5: Revenue of harvest
+* [x] Include cost function for module 6: Penalties in case of exceeding legal limits
+* [ ] Include cost function for module 7: Cost of greenhouse gas emissions
+* [ ] Include realistic cultivation response curves from module 5
+* [ ] Add other nutrients than NPK to module 5
+* [x] Add function to update fertilizer properties
 * [ ] Add custom fertilizers
 * [ ] Add custom cultivations
 * [ ] Add details of the optimization to the result

data-raw/fertilizers.R

@@ -13,5 +13,28 @@ token <- ''
 fertilizers <- fread(paste0('https://raw.githubusercontent.com/AgroCares/pandex/main/data-raw/b_fp/b_fp_srm.csv?token=', token))
 setnames(fertilizers, colnames(fertilizers), tolower(colnames(fertilizers)))
 
+
+# Assign an id ------------------------------------------------------------
+fertilizers[, p_id := 1:.N]
+
+
+# Assign derivative parameters --------------------------------------------
+fertilizers[, p_stored := 0]
+fertilizers[, p_price := 1]
+
+fertilizers[, p_type_manure := fifelse(p_type_manure, 1, 0)]
+fertilizers[, p_p_wcl := 1]
+fertilizers[p_type_compost == TRUE, p_p_wcl := 0.25]
+fertilizers[p_name_nl == 'Champost', p_p_wcl := 0.75]
+fertilizers[p_name_nl == 'Rundvee vaste mest', p_p_wcl := 0.75]
+
+fertilizers[, p_storage_cost := 10000]
+fertilizers[, p_storage_capacity := 1000000]
+fertilizers[, p_storage_available := 0]
+fertilizers[p_type_artificial == TRUE, p_storage_capacity := 1000]
+fertilizers[p_type_artificial == TRUE, p_storage_cost := 100]
+fertilizers[p_type_artificial == TRUE, p_storage_available := 1]
+
 # Export table ------------------------------------------------------------
 usethis::use_data(fertilizers, overwrite = TRUE, version = 3, compress = 'xz')
+