From 7b9d9b5fc930b43da9ec0423a9f232d68932835e Mon Sep 17 00:00:00 2001
From: Andreas Gattringer <andreas.gattringer@univie.ac.at>
Date: Fri, 23 Jun 2023 23:30:35 +0200
Subject: [PATCH] butterflies: vital rates changes

---
 .../butterflies/butterflies_generations.c     |  8 ++++-
 src/modules/butterflies/butterflies_main.c    | 12 ++++---
 src/modules/butterflies/butterflies_main.h    |  2 +-
 .../butterflies/butterflies_overlays.c        |  6 ++++
 .../butterflies/butterflies_populations.c     | 35 +++++++++----------
 .../butterflies/butterflies_populations.h     |  2 +-
 .../butterflies/butterflies_vital_rates.c     | 23 ++++--------
 7 files changed, 46 insertions(+), 42 deletions(-)

diff --git a/src/modules/butterflies/butterflies_generations.c b/src/modules/butterflies/butterflies_generations.c
index c2267e7..8252f0c 100644
--- a/src/modules/butterflies/butterflies_generations.c
+++ b/src/modules/butterflies/butterflies_generations.c
@@ -21,12 +21,18 @@ void bf_area_generation_update(struct cats_grid *grid, struct cats_thread_info *
         const cats_dt_coord start_col = ts->area.start_col;
         const cats_dt_coord end_col = ts->area.end_col;
 
+        cats_dt_rates ot = grid->param.OT;
+        bool hybrid = grid->param.parametrization == PARAM_HYBRID;
 
         for (cats_dt_coord row = start_row; row < end_row; row++) {
                 for (cats_dt_coord col = start_col; col < end_col; col++) {
+                        bool suitability_to_low = false;
+
+                        if (hybrid && get_suitability(grid, row, col) < ot) suitability_to_low = true;
 
                         if (cell_excluded_by_overlay(conf, row, col)
-                            ||  ! (data->info_layer[row][col] & BF_CELL_VALID_DISPERSAL_TARGET)) {
+                            ||  ! (data->info_layer[row][col] & BF_CELL_VALID_DISPERSAL_TARGET)
+                            ||  suitability_to_low == true) {
                                 data->eggs[row][col] = 0.0f;
                                 data->generations[row][col] = 0.0f;
                                 set_population_ignore_cc(grid, row, col, 0);
diff --git a/src/modules/butterflies/butterflies_main.c b/src/modules/butterflies/butterflies_main.c
index e9650a4..3c6583a 100644
--- a/src/modules/butterflies/butterflies_main.c
+++ b/src/modules/butterflies/butterflies_main.c
@@ -24,17 +24,19 @@ double *butterflies_matrix(struct cats_configuration *conf, struct lambda_parame
 
         struct conf_data_butterflies *module_conf = CATS_MODULE_DATA;
         cats_dt_rates eggs_per_female = calculate_rate_for_matrix(&module_conf->eggs_per_female, l_param, print_rate);
-        cats_dt_rates adults_per_female = calculate_rate_for_matrix(&module_conf->adults_per_female, l_param, print_rate);
+        cats_dt_rates reproduction_rate = calculate_rate_for_matrix(&module_conf->reproduction_rate, l_param, print_rate);
         cats_dt_rates female_fraction = module_conf->female_fraction;
         cats_dt_rates stationary = module_conf->probability_to_stay;
         cats_dt_rates mobile = 1.0 - stationary;
         cats_dt_rates egg_fraction_source = module_conf->egg_fraction_source;
-        cats_dt_rates max_eggs = module_conf->eggs_per_female.max_rate;
 
-        cats_dt_rates eggs_to_adults_rate = bf_egg_to_adult_survival_rate(adults_per_female, max_eggs);
         cats_dt_rates local_eggs =  (stationary + mobile * egg_fraction_source ) * eggs_per_female;
-        // female -> female;
+        // female -> female
+        // to achieve the target reproduction rate, the number of eggs per female laid in the cell
+        // that survive and become adult has to be the reproduction rate divided by the female fraction divided by the number of eggs
+        cats_dt_rates eggs_to_adults_rate = bf_egg_to_adult_survival_rate(reproduction_rate, local_eggs, female_fraction);
         cats_dt_rates result =  local_eggs * eggs_to_adults_rate * female_fraction;
+        log_message(LOG_INFO, "scale factor %Lf: eggs_to_adults_rate %Lf", l_param->grid->param.scale_factor, eggs_to_adults_rate);
 
         double *matrix = calloc_or_die(1, sizeof(double));
 
@@ -95,7 +97,7 @@ void load_butterflies_species_params(struct cats_configuration *conf, struct cat
 
         load_conf_vital_rate(&data->eggs_per_female, conf, ini, section_name, param);
         load_conf_vital_rate(&data->butterfly_generations, conf, ini, section_name, param);
-        load_conf_vital_rate(&data->adults_per_female, conf, ini, section_name, param);
+        load_conf_vital_rate(&data->reproduction_rate, conf, ini, section_name, param);
 
 
         //load_conf_vital_rate(&data->butterfly_egg_to_adult_survival, conf, ini, section_name, param);
diff --git a/src/modules/butterflies/butterflies_main.h b/src/modules/butterflies/butterflies_main.h
index 9da1de6..b122556 100644
--- a/src/modules/butterflies/butterflies_main.h
+++ b/src/modules/butterflies/butterflies_main.h
@@ -47,7 +47,7 @@ struct conf_data_butterflies {
     struct cats_vital_rate eggs_per_female;
     struct cats_vital_rate butterfly_egg_to_adult_survival;
     struct cats_vital_rate butterfly_generations;
-    struct cats_vital_rate adults_per_female;
+    struct cats_vital_rate reproduction_rate;
     cats_dt_rates female_fraction;
 
     bool actions_added;
diff --git a/src/modules/butterflies/butterflies_overlays.c b/src/modules/butterflies/butterflies_overlays.c
index 9a937c2..4e36332 100644
--- a/src/modules/butterflies/butterflies_overlays.c
+++ b/src/modules/butterflies/butterflies_overlays.c
@@ -6,6 +6,7 @@
 #include "butterflies_overlays.h"
 #include "butterflies_main.h"
 #include "inline_overlays.h"
+#include "inline_population.h"
 
 
 enum action_status bf_grid_overlay_update(const struct cats_configuration *conf, struct cats_grid *grid)
@@ -32,6 +33,8 @@ enum action_status bf_grid_overlay_update(const struct cats_configuration *conf,
 
                         if (cell_excluded_by_overlay(conf, row, col)) {
                                 data->info_layer[row][col] |= BF_CELL_EXCLUDED;
+                                data->eggs[row][col] = 0;
+                                set_population_ignore_cc(grid, row, col, 0);
                                 cells_excluded += 1;
                                 continue;
                         }
@@ -54,6 +57,9 @@ enum action_status bf_grid_overlay_update(const struct cats_configuration *conf,
                                 data->info_layer[row][col] |= BF_CELL_VALID_DISPERSAL_TARGET;
                                 cells_habitat_and_resource_ok += 1;
 
+                        } else {
+                                data->eggs[row][col] = 0;
+                                set_population_ignore_cc(grid, row, col, 0);
                         }
                 }
         }
diff --git a/src/modules/butterflies/butterflies_populations.c b/src/modules/butterflies/butterflies_populations.c
index a5e9755..5d40f45 100644
--- a/src/modules/butterflies/butterflies_populations.c
+++ b/src/modules/butterflies/butterflies_populations.c
@@ -17,12 +17,9 @@
 #include "populations/population.h"
 
 
-cats_dt_rates bf_egg_to_adult_survival_rate(cats_dt_rates adults_per_female, cats_dt_rates max_eggs)
+cats_dt_rates bf_egg_to_adult_survival_rate(cats_dt_rates reproduction_rate, cats_dt_rates eggs, cats_dt_rates female_fraction)
 {
-
-
-        return adults_per_female/max_eggs;
-
+        return (reproduction_rate / female_fraction) / eggs;
 }
 
 void bf_cell_maturation(struct cats_grid *grid, struct cats_thread_info *ts, cats_dt_coord row, cats_dt_coord col,
@@ -62,7 +59,7 @@ void bf_cell_maturation(struct cats_grid *grid, struct cats_thread_info *ts, cat
         // first the all cells which have 5 (> 4) generations are modelled
         // then all the cells with have at least 4 (> 3) generations are modelled
         // ...
-        // at last all cells with have at least 1 generations are modelled
+        // at last all cells with have at least 1 generation are modelled
 
         float this_generation_fraction = 1.0f;
         int32_t local_max_generation = (int32_t) ceilf(data->generations[row][col]);
@@ -89,25 +86,27 @@ void bf_cell_maturation(struct cats_grid *grid, struct cats_thread_info *ts, cat
 
         // not capped, we can have more adults than CC
         struct conf_data_butterflies *module_conf = CATS_MODULE_DATA;
-        cats_dt_rates adults_per_female = calculate_rate(&module_conf->adults_per_female, NAN, &grid->param,
+        cats_dt_environment suit = get_suitability(grid, row, col);
+        cats_dt_rates reproduction_rate = calculate_rate(&module_conf->reproduction_rate, NAN, &grid->param,
                                                          grid, row, col, NULL);
+        if (suit < grid->param.OT && reproduction_rate > 0) {
+                log_message(LOG_ERROR, "Suitability %f under OT %Lf, but adults per female = %Lf", suit, grid->param.OT, reproduction_rate);
+                exit_cats(EXIT_FAILURE);
+        }
 
-        cats_dt_rates max_eggs = module_conf->eggs_per_female.max_rate;
-        cats_dt_rates survival = bf_egg_to_adult_survival_rate(adults_per_female, max_eggs);
+        if (reproduction_rate == 0) {
+                data->eggs[row][col] = 0;
+                set_population_ignore_cc(grid, row, col, 0);
+                return;
+        };
+
+        cats_dt_rates survival = bf_egg_to_adult_survival_rate(reproduction_rate, eggs, module_conf->female_fraction);
         cats_dt_population adults = poisson(ts->rng, eggs * survival);
         assert(adults >= 0);
-        //cats_dt_environment suit = get_suitability(grid, row, col);
-        //printf("ABCD,%d,%d,%0.3f\n", (int) orig_eggs, adults, suit);
+
 
         set_population(grid, row, col, adults);
-/*
-        if (adults > max_cc * 10) {
 
-                log_message(LOG_ERROR, "row %d col %d: number of adults %d exceeds %d time smaximum carrying capacity %d (maturation %Lf, eggs %f, max adults per f %Lf), suit %f, OT %Lf, survival %Lf", row, col, adults, 10, max_cc,
-                            survival, eggs, module_conf->adults_per_female.max_rate, suit, grid->param.OT, survival);
-                exit_cats(EXIT_FAILURE);
-        }
-        */
 }
 
 
diff --git a/src/modules/butterflies/butterflies_populations.h b/src/modules/butterflies/butterflies_populations.h
index eed6c10..f8f164b 100644
--- a/src/modules/butterflies/butterflies_populations.h
+++ b/src/modules/butterflies/butterflies_populations.h
@@ -7,6 +7,6 @@
 #include "data/cats_grid.h"
 void bf_cell_maturation(struct cats_grid *grid, struct cats_thread_info *ts, cats_dt_coord row, cats_dt_coord col, bool check_exclusion);
 void bf_area_kill_adults(struct cats_grid *grid, struct cats_thread_info *ts);
-cats_dt_rates bf_egg_to_adult_survival_rate(cats_dt_rates adults_per_female, cats_dt_rates max_eggs);
+cats_dt_rates bf_egg_to_adult_survival_rate(cats_dt_rates reproduction_rate, cats_dt_rates eggs, cats_dt_rates female_ratio);
 
 #endif //CATS_BUTTERFLIES_POPULATIONS_H
diff --git a/src/modules/butterflies/butterflies_vital_rates.c b/src/modules/butterflies/butterflies_vital_rates.c
index 307a7e8..de12309 100644
--- a/src/modules/butterflies/butterflies_vital_rates.c
+++ b/src/modules/butterflies/butterflies_vital_rates.c
@@ -2,7 +2,7 @@
 
 void bf_add_vital_rates(struct cats_configuration *conf, struct conf_data_butterflies *data)
 {
-        // generations
+        // number of generations
         register_module_vital_rate(conf, &data->butterfly_generations,"butterflies generations");
         set_vital_rate_name(&data->butterfly_generations, "butterflies generations");
         set_vital_rate_suitability_cutoff_hint(&data->butterfly_generations, HYBRID_SUIT_TS_ZT);
@@ -13,21 +13,12 @@ void bf_add_vital_rates(struct cats_configuration *conf, struct conf_data_butter
         // adult to eggs
         register_module_vital_rate(conf, &data->eggs_per_female, "butterflies eggs per female");
         set_vital_rate_name(&data->eggs_per_female, "butterflies eggs per female");
-        set_vital_rate_suitability_cutoff_hint(&data->eggs_per_female, HYBRID_SUIT_TS_ZT);
-        set_vital_rate_link_hybrid_function(&data->eggs_per_female, conf, LINK_SUITABILITY_SIGMOID);
+        set_vital_rate_suitability_cutoff_hint(&data->eggs_per_female, HYBRID_SUIT_TS_OT);
+        set_vital_rate_link_hybrid_function(&data->eggs_per_female, conf, LINK_CONSTANT);
 
         // eggs to adults
-        register_module_vital_rate(conf, &data->adults_per_female, "butterflies adults per female");
-        set_vital_rate_name(&data->adults_per_female, "butterflies adults per female");
-        set_vital_rate_suitability_cutoff_hint(&data->adults_per_female, HYBRID_SUIT_TS_OT);
-        set_vital_rate_link_hybrid_function(&data->adults_per_female, conf, LINK_SUITABILITY_SIGMOID);
-
-        // egg to adult survival
-        /*
-        register_module_vital_rate(conf, &data->butterfly_egg_to_adult_survival,"butterflies egg to adult survival rate");
-        set_vital_rate_name(&data->butterfly_egg_to_adult_survival, "butterflies eggs to adult survival rate");
-        set_vital_rate_suitability_cutoff_hint(&data->butterfly_egg_to_adult_survival, HYBRID_SUIT_TS_ZT);
-        set_vital_rate_link_hybrid_function(&data->butterfly_egg_to_adult_survival, conf, LINK_SUITABILITY_SIGMOID);
-        set_vital_density(&data->butterfly_egg_to_adult_survival, NO_DENSITY_DEP);
-         */
+        register_module_vital_rate(conf, &data->reproduction_rate, "butterflies reproduction rate");
+        set_vital_rate_name(&data->reproduction_rate, "butterflies reproduction rate");
+        set_vital_rate_link_hybrid_function(&data->reproduction_rate, conf, LINK_SUITABILITY_SIGMOID);
+        set_vital_rate_suitability_cutoff_hint(&data->reproduction_rate, HYBRID_SUIT_TS_OT);
 }
-- 
GitLab