From a7fc88a36dc27262b67eee97a878337416065896 Mon Sep 17 00:00:00 2001 From: Rafael Sachetto Date: Mon, 19 Aug 2024 11:29:44 -0300 Subject: [PATCH] Adding simple stimulus support to the eikonal solver --- example_configs/rabbit_mesh_eikonal.ini | 16 +++++++ src/main_eikonal.c | 60 +++++++++++++++++++++---- src/stimuli_library/stimuli.c | 42 +++++++++++++++++ 3 files changed, 109 insertions(+), 9 deletions(-) diff --git a/example_configs/rabbit_mesh_eikonal.ini b/example_configs/rabbit_mesh_eikonal.ini index 6df592a2..45454435 100644 --- a/example_configs/rabbit_mesh_eikonal.ini +++ b/example_configs/rabbit_mesh_eikonal.ini @@ -11,6 +11,22 @@ main_function=initialize_grid_with_rabbit_mesh ;These can be optional depending on the domain main_function mesh_file=meshes/rabheart.alg +[stim_point] +current = 1 +center_x=31875 +center_y=10875 +center_z=16375 +start = 0.0 +duration = 2.0 +main_function= stim_if_point_equal + +;[stim_mouse] +;x_limit=500.0 +;start = 0.0 +;duration = 2.0 +;current = -50.0 +;main_function= stim_if_x_less_than + [save_result] ;/////mandatory///////// output_dir=./outputs/rabbit_eikonal diff --git a/src/main_eikonal.c b/src/main_eikonal.c index 61d7418c..26ec13f2 100644 --- a/src/main_eikonal.c +++ b/src/main_eikonal.c @@ -5,6 +5,7 @@ #include "config/config_parser.h" #include "config/domain_config.h" #include "eikonal/eikonal_solver.h" +#include "config/stim_config.h" #include "logger/logger.h" #include "utils/file_utils.h" #include @@ -19,7 +20,7 @@ static int compare_coordinates(const void *a, const void *b) { } if (coord1->center.x != coord2->center.x) { - return (int) (coord1->center.x - coord2->center.x); + return (int) (coord1->center.x - coord2->center.x); } @@ -32,6 +33,9 @@ int main(int argc, char *argv[]) { struct eikonal_options *eikonal_options = new_eikonal_options(); parse_eikonal_options(argc, argv, eikonal_options); struct grid *grid = new_grid(); + struct string_voidp_hash_entry *stimuli_configs = eikonal_options->stim_configs; + + struct eikonal_solver *solver = new_eikonal_solver(false); if(ini_parse(eikonal_options->config_file, parse_eikonal_config_file, eikonal_options) < 0) { log_error_and_exit("Error: Can't load the config file %s\n", eikonal_options->config_file); @@ -59,6 +63,52 @@ int main(int argc, char *argv[]) { } + if(stimuli_configs) { + + size_t n = shlen(stimuli_configs); + struct time_info time_info = {0.0, 0.0, 0.0, 0}; + + if(n > 0) { + STIM_CONFIG_HASH_FOR_INIT_FUNCTIONS(stimuli_configs); + } + + for(int i = 0; i < n; i++) { + + struct string_voidp_hash_entry e = stimuli_configs[i]; + log_info("Stimulus name: %s\n", e.key); + struct config *tmp = (struct config *)e.value; + print_stim_config_values(tmp); + log_msg(LOG_LINE_SEPARATOR); + } + + set_spatial_stim(&time_info, stimuli_configs, grid, false); + + struct config *tmp = NULL; + + uint32_t n_active = grid->num_active_cells; + struct cell_node **ac = grid->active_cells; + + for(size_t k = 0; k < n; k++) { + + tmp = (struct config *)stimuli_configs[k].value; + + for(uint32_t i = 0; i < n_active; i++) { + real value = ((real *)(tmp->persistent_data))[i]; + if(value != 0.0) { + size_t *initial_seed = calloc(sizeof(size_t), 3); + double new_pos_x = (ac[i]->center.x-ac[i]->discretization.x/2.0)/ac[i]->discretization.x; + double new_pos_y = (ac[i]->center.y-ac[i]->discretization.y/2.0)/ac[i]->discretization.y; + double new_pos_z = (ac[i]->center.z-ac[i]->discretization.z/2.0)/ac[i]->discretization.z; + + initial_seed[0] = new_pos_x; + initial_seed[1] = new_pos_y; + initial_seed[2] = new_pos_z; + arrput(solver->seeds, initial_seed); + } + } + } + } + //Translating to 0 for(int i = 0 ; i < grid->num_active_cells; i++) { double new_pos_x = (grid->active_cells[i]->center.x-grid->active_cells[i]->discretization.x/2.0)/grid->active_cells[i]->discretization.x; @@ -76,7 +126,6 @@ int main(int argc, char *argv[]) { size_t itersPerBlock = 10, type = 1; - struct eikonal_solver *solver = new_eikonal_solver(false); solver->width = (int) (grid->cube_side_length.x/grid->active_cells[0]->discretization.x); solver->height = (int) (grid->cube_side_length.y/grid->active_cells[0]->discretization.y); solver->depth = (int) (grid->cube_side_length.z/grid->active_cells[0]->discretization.z); @@ -84,12 +133,6 @@ int main(int argc, char *argv[]) { solver->solver_type = type; solver->iters_per_block = itersPerBlock; - size_t *initial_seed = calloc(sizeof(size_t), 3); - initial_seed[0] = 128; - initial_seed[1] = 48; - initial_seed[2] = 63; - - arrput(solver->seeds, initial_seed); solver->active_cells = grid->active_cells; solver->num_active_cells = grid->num_active_cells; @@ -115,7 +158,6 @@ int main(int argc, char *argv[]) { } } - free(initial_seed); free_eikonal_solver(solver); free_eikonal_options(eikonal_options); diff --git a/src/stimuli_library/stimuli.c b/src/stimuli_library/stimuli.c index 75797fbd..8adadcae 100644 --- a/src/stimuli_library/stimuli.c +++ b/src/stimuli_library/stimuli.c @@ -11,6 +11,48 @@ #include "../config_helpers/config_helpers.h" #include "../utils/utils.h" +SET_SPATIAL_STIM(stim_if_point_equal) { + + uint32_t n_active = the_grid->num_active_cells; + struct cell_node **ac = the_grid->active_cells; + + if(is_purkinje) { + n_active = the_grid->purkinje->num_active_purkinje_cells; + ac = the_grid->purkinje->purkinje_cells; + } + + ALLOCATE_STIMS(); + + real stim_current = 0.0; + GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real, stim_current, config, "current"); + + bool stim; + real stim_value; + + real_cpu center_x = 0.0; + real_cpu center_y = 0.0; + real_cpu center_z = 0.0; + GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, center_x, config, "center_x"); + GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, center_y, config, "center_y"); + GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, center_z, config, "center_z"); + + uint32_t i; + + OMP(parallel for private(stim, stim_value)) + for(i = 0; i < n_active; i++) { + stim = ac[i]->center.x == center_x && ac[i]->center.y == center_y && ac[i]->center.z == center_z; + + if(stim) { + stim_value = stim_current; + } else { + stim_value = 0.0; + } + + SET_STIM_VALUE(i, stim_value); + } + +} + SET_SPATIAL_STIM(stim_if_x_less_than) { uint32_t n_active = the_grid->num_active_cells;