Skip to content
Snippets Groups Projects
Commit 0a2d167e authored by Andreas Gattringer's avatar Andreas Gattringer
Browse files

updated module calculation of lambda/scale factor

parent 3b2c0ead
Branches
Tags
No related merge requests found
......@@ -42,34 +42,47 @@
#include "lambda/leslie_matrix.h"
#include "lambda/lambda.h"
struct cats_module *get_demographic_module_for_matrix(struct cats_configuration *conf, int32_t species_id)
{
struct cats_module *module = NULL;
struct cats_species_param *param = &conf->param[species_id];
if (param->demographic_module == NULL) return module;
int32_t module_id = param->demographic_module->id;
if (conf->modules.module[module_id].flags & MODULE_ALTERNATE_DEMOGRAPHIC) return &conf->modules.module[module_id];
return NULL;
}
bool have_leslie_matrix_function_for_module(struct cats_module *module)
{
if (module == NULL) return true;
cats_leslie_matrix_func func = module->create_leslie_matrix;
return func != NULL;
}
cats_dt_rates calculate_scale_factor_from_lambda(struct cats_configuration *conf, int32_t species_id)
{
assert(species_id < conf->param_count);
struct cats_species_param *param = &conf->param[species_id];
bool ok = false;
cats_leslie_matrix_func func = NULL;
// FIXME LOOKS LIKE A MESS
if (param->custom_vital_rates == false || param->custom_vital_ages == false) ok = true;
if (param->custom_vital_rates || param->custom_vital_ages) {
if (param->demographic_module != NULL) {
int32_t id = param->demographic_module->id;
if (conf->modules.module[id].flags & MODULE_ALTERNATE_DEMOGRAPHIC) {
func = conf->modules.module[id].create_leslie_matrix;
if (func != NULL) {
ok = true;
}
}
}
}
struct cats_module *module = get_demographic_module_for_matrix(conf, species_id);
bool have_func = have_leslie_matrix_function_for_module(module);
log_message(LOG_IMPORTANT, "Demographic module %p", module);
if (module != NULL && have_func == false) {
if (!ok) {
log_message(LOG_WARNING, "unable to calculate scale factor for %s: custom vital rates or vital ages",
param->species_name);
log_message(LOG_WARNING, "Using scale factor 0.5 for %s", param->species_name);
return 0.5;
} else if (module != NULL) {
log_message(LOG_IMPORTANT, "Grid %d uses module %s for demography: %p", species_id,
module->canonical_name, module->create_leslie_matrix);
}
......@@ -104,8 +117,9 @@ cats_dt_rates calculate_scale_factor_from_lambda(struct cats_configuration *conf
conf->param[species_id].scale_factor = s;
conf->direct_scale_factor = s;
cats_dt_rates old_s = s;
lambda = calculate_lambda(conf, &l_param, true, func);
log_message(LOG_INFO, "iteration %d: scale factor %0.25Lf resulted in lambda %0.25Lf", iterations, s, lambda);
lambda = calculate_lambda(conf, &l_param, true, module);
log_message(LOG_INFO, "iteration %d: scale factor %0.25Lf resulted in lambda %0.25Lf", iterations, s,
lambda);
if (fabsl(lambda - 1.0) < threshold) {
conf->lambda_at_OT = lambda;
......@@ -164,7 +178,7 @@ cats_dt_rates calculate_scale_factor_from_lambda(struct cats_configuration *conf
iterations);
conf->param[species_id].scale_factor = s;
conf->direct_scale_factor = s;
calculate_lambda(conf, &l_param, false, func);
calculate_lambda(conf, &l_param, false, module);
return s;
......
......@@ -78,13 +78,15 @@ lookup_lambda(struct cats_configuration *conf, struct cats_grid *grid, cats_dt_p
l_param.suitability = NAN;
break;
}
/*
cats_leslie_matrix_func func = NULL;
if (grid->param.demographic_module && grid->param.demographic_module->create_leslie_matrix) {
func = grid->param.demographic_module->create_leslie_matrix;
}
*/
// we have density dependence
if (N >= 0 || grid->param.parametrization != PARAM_HYBRID) {
return calculate_lambda(conf, &l_param, true, func);
return calculate_lambda(conf, &l_param, true, NULL);
}
assert(grid->param.parametrization == PARAM_HYBRID);
......@@ -100,7 +102,7 @@ lookup_lambda(struct cats_configuration *conf, struct cats_grid *grid, cats_dt_p
return param->lambda_cache[index];
} else {
l_param.N = 0;
cats_dt_rates lambda = calculate_lambda(conf, &l_param, true, func);
cats_dt_rates lambda = calculate_lambda(conf, &l_param, true, NULL);
param->lambda_cache[index] = lambda;
return lambda;
}
......
......@@ -39,18 +39,29 @@ void debug_lambda_problem(struct cats_configuration *conf, gsl_complex largest_e
cats_dt_rates calculate_lambda(struct cats_configuration *conf, struct lambda_parameters *l_param, bool silent,
cats_leslie_matrix_func leslie_matrix_func)
struct cats_module *module)
{
int32_t N = 0;
bool debug = conf->command_line_options.debug_flags & DEBUG_LAMBDA;
const int32_t species_id = l_param->species_id;
double *matrix = NULL;
if (leslie_matrix_func == NULL) {
int32_t N = 0;
if (module == NULL) {
matrix = create_leslie_matrix(conf, l_param, silent, &N);
} else if (module->create_leslie_matrix != NULL) {
matrix = module->create_leslie_matrix(conf, l_param, silent, &N);
} else {
matrix = leslie_matrix_func(conf, l_param, silent, &N);
log_message(LOG_ERROR, "no function to create leslie matrix specified");
exit_cats(EXIT_FAILURE);
}
bool debug = conf->command_line_options.debug_flags & DEBUG_LAMBDA;
const int32_t species_id = l_param->species_id;
//
double *matrix_copy = NULL;
......@@ -76,18 +87,22 @@ cats_dt_rates calculate_lambda(struct cats_configuration *conf, struct lambda_pa
print_matrix(matrix_copy, N, N);
printf("Schur form T:\n");
print_matrix(matrix, N, N);
if (module == NULL) {
debug_lambda_1_found(conf, species_id, largest_eigen_value, matrix_copy, N, l_param);
printf("============================\n");
}
}
if (silent == false && l_param->calculate_scale == true) {
if (fabs(GSL_REAL(largest_eigen_value) - 1.0) < LAMBDA_EPS) {
printf("Final matrix\n");
print_matrix(matrix_copy, N, N);
if (module == NULL) {
debug_lambda_1_found(conf, species_id, largest_eigen_value, matrix_copy, N, l_param);
}
}
}
// FREE EVERYTHING
gsl_matrix_complex_free(eigen_vectors);
......
......@@ -28,6 +28,6 @@
#include "leslie_matrix.h"
cats_dt_rates calculate_lambda(struct cats_configuration *conf, struct lambda_parameters *l_param, bool silent,
cats_leslie_matrix_func leslie_matrix_func);
struct cats_module *module);
#endif //CATS_LAMBDA_H
......@@ -36,6 +36,37 @@
#include "inline_vital_rates.h"
cats_dt_rates
calculate_rate_for_matrix(struct cats_vital_rate *rate,
const struct lambda_parameters *l_param,
bool print)
{
assert(l_param != NULL);
assert(l_param->param != NULL);
const struct cats_species_param *param = l_param->param; //&conf->param[grid_id];
const cats_dt_rates N = l_param->N;
const cats_dt_rates K = l_param->K;
struct link_override_parameters override = {-1};
if (param->parametrization == PARAM_HYBRID) {
override.override_suitability = l_param->suitability;
}
override.override_population = N;
override.override_carrying_capacity = l_param->K;
cats_dt_rates value = calculate_rate(rate, N, param, l_param->grid, l_param->row, l_param->col, &override);
if (print) {
log_message(LOG_INFO, "%s: %Lf (maximum %Lf) for population %Lf/%Lf",
rate->name, value, rate->max_rate, N, K);
}
return value;
}
cats_dt_rates
get_rate_for_matrix(enum cats_vital_rate_id rate_type,
const struct lambda_parameters *l_param,
......
......@@ -43,6 +43,10 @@ cats_dt_rates
get_rate_for_matrix(enum cats_vital_rate_id rate_type, const struct lambda_parameters *l_param,
bool print);
cats_dt_rates
calculate_rate_for_matrix(struct cats_vital_rate *rate,
const struct lambda_parameters *l_param,
bool print);
double *
create_leslie_matrix(struct cats_configuration *conf, struct lambda_parameters *l_param, bool silent, int32_t *N_out);
......
......@@ -148,8 +148,8 @@ void register_load_species_param_config_func(struct cats_configuration *conf, ca
void register_create_leslie_matrix_func(struct cats_configuration *conf, cats_leslie_matrix_func func)
{
conf->param[0].module_data[CATS_MODULE_ID].create_leslie_matrix_function = func;
conf->modules.module[CATS_MODULE_ID].create_leslie_matrix = func;
printf("LESLIE MATRIX FUNC: %p\n", conf->modules.module[CATS_MODULE_ID].create_leslie_matrix);
}
......
......@@ -28,7 +28,6 @@
#include "modules.h"
#include "cats_strings/cats_strings.h"
void init_cats_module(struct cats_module *module, int32_t id)
{
module->name = NULL;
......@@ -38,6 +37,8 @@ void init_cats_module(struct cats_module *module, int32_t id)
module->id = id;
module->create_leslie_matrix = NULL;
module->stats_header = NULL;
}
......@@ -59,7 +60,6 @@ void init_module_species_data(struct module_species_data *m)
m->grid_init_function = NULL;
m->grid_cleanup_function = NULL;
m->load_species_param_function = NULL;
m->create_leslie_matrix_function = NULL;
m->grid_stat_function = NULL;
m->module = NULL;
m->vital_rate_count = 0;
......
......@@ -23,12 +23,11 @@
#ifndef CATS_MODULES_H
#define CATS_MODULES_H
#include <stdint.h>
#include <stdbool.h>
#include "defaults.h"
#include "data/cats_datatypes.h"
#include <gsl/gsl_eigen.h>
struct cats_ini;
struct cats_configuration;
struct lambda_parameters;
......@@ -53,13 +52,14 @@ typedef struct string_array *(*cats_grid_stat_function)(struct cats_configuratio
typedef void (*cats_cell_function)(struct cats_configuration *conf, struct cats_grid *grid, cats_dt_coord row,
cats_dt_coord col, cats_dt_population K);
typedef void (*eigen_system_print_func)(gsl_complex eigen_value, gsl_vector_complex_view *eigen_vector, int32_t N, struct cats_species_param *param);
struct module_species_data {
struct cats_module *module;
cats_grid_function grid_init_function;
cats_grid_function grid_cleanup_function;
cats_load_species_param_function load_species_param_function;
cats_leslie_matrix_func create_leslie_matrix_function;
cats_grid_stat_function grid_stat_function;
cats_cell_function cell_destroyed_action;
cats_cell_function cell_carrying_capacity_action;
......@@ -68,7 +68,6 @@ struct module_species_data {
char **vital_rate_name;
bool *vital_rate_required;
int32_t vital_rate_count;
void *module_data;
};
......@@ -79,13 +78,16 @@ enum cats_module_flags {
};
struct cats_module {
char *name;
char *canonical_name;
char *filename;
void *data;
enum cats_module_flags flags;
void *create_leslie_matrix;
cats_leslie_matrix_func create_leslie_matrix;
int32_t id;
struct string_array *stats_header;
};
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment