Merge pull request #68 from bergolho/master

rsachetto · web-flow · commit f461dee1b376 · 2024-04-24T10:05:46.000-03:00
New conic domain function and more general format for the benchmark domain function
diff --git a/src/domains_library/domain.c b/src/domains_library/domain.c
@@ -127,6 +127,10 @@ SET_SPATIAL_DOMAIN(initialize_grid_with_rabbit_mesh) {
     return num_loaded > 0;
 }
 
+/**
+ * Sets the current domain as a domain described in the N-version benchmark
+ * (http://rsta.royalsocietypublishing.org/content/369/1954/4331)
+ */
 SET_SPATIAL_DOMAIN(initialize_grid_with_benchmark_mesh) {
 
     real_cpu side_length;
@@ -137,6 +141,15 @@ SET_SPATIAL_DOMAIN(initialize_grid_with_benchmark_mesh) {
     real_cpu max_h = start_h;
     GET_PARAMETER_NUMERIC_VALUE_OR_USE_DEFAULT(real_cpu, max_h, config, "maximum_discretization");
 
+    real_cpu side_length_x = 20000.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, side_length_x, config, "side_length_x");
+
+    real_cpu side_length_y = 7000.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, side_length_y, config, "side_length_y");
+
+    real_cpu side_length_z = 3000.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, side_length_z, config, "side_length_z");
+
     log_info("Loading N-Version benchmark mesh using dx %lf um, dy %lf um, dz %lf um\n", start_h, start_h, start_h);
 
     side_length = start_h;
@@ -150,7 +163,7 @@ SET_SPATIAL_DOMAIN(initialize_grid_with_benchmark_mesh) {
     int num_steps = get_num_refinement_steps_to_discretization(side_length, start_h);
 
     refine_grid(the_grid, num_steps);
-    set_benchmark_domain(the_grid);
+    set_benchmark_domain(the_grid, side_length_x, side_length_y, side_length_z);
 
     log_info("Cleaning grid\n");
     int i;
@@ -401,5 +414,52 @@ SET_SPATIAL_DOMAIN(initialize_grid_with_square_mesh_and_source_sink_fibrotic_reg
 
     set_plain_fibrosis_source_sink_region(the_grid, phi, seed, min_x, max_x, min_y, max_y, min_z, max_z, source_sink_min_x, source_sink_max_x, side_length);
 
+    return 1;
+}
+
+SET_SPATIAL_DOMAIN(initialize_grid_with_cuboid_and_sphere_fibrotic_mesh_with_conic_path){
+    real_cpu side_length_x = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, side_length_x, config, "side_length_x");
+
+    real_cpu side_length_y = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, side_length_y, config, "side_length_y");
+
+    real_cpu side_length_z = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, side_length_z, config, "side_length_z");
+
+    real_cpu start_dx = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, start_dx, config, "start_dx");
+
+    real_cpu start_dy = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, start_dy, config, "start_dy");
+
+    real_cpu start_dz = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, start_dz, config, "start_dz");
+
+    real_cpu phi = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, phi, config, "phi");
+
+    real_cpu plain_center = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, plain_center, config, "plain_center");
+
+    real_cpu sphere_radius = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, sphere_radius, config, "sphere_radius");
+
+    real_cpu border_zone_size = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, border_zone_size, config, "border_zone_size");
+
+    real_cpu border_zone_radius = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, border_zone_radius, config, "border_zone_radius");
+
+    real_cpu conic_slope = 0.0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_REPORT_ERROR(real_cpu, conic_slope, config, "conic_slope");
+
+    unsigned seed = 0;
+    GET_PARAMETER_NUMERIC_VALUE_OR_USE_DEFAULT(unsigned, seed, config, "seed");
+
+    //set_square_mesh(config, the_grid);
+    set_cuboid_domain_mesh(the_grid, start_dx, start_dy, start_dz, side_length_x, side_length_y, side_length_z);
+    set_cuboid_sphere_fibrosis_with_conic_path(the_grid, phi, plain_center, sphere_radius, border_zone_size, border_zone_radius, seed, conic_slope);
+
     return 1;
 }
diff --git a/src/domains_library/domain_helpers.c b/src/domains_library/domain_helpers.c
@@ -445,13 +445,13 @@ int calculate_cuboid_side_lengths(real_cpu start_dx, real_cpu start_dy, real_cpu
  * (http://rsta.royalsocietypublishing.org/content/369/1954/4331)
  *
  */
-void set_benchmark_domain(struct grid *the_grid) {
+void set_benchmark_domain(struct grid *the_grid, real_cpu sx, real_cpu sy, real_cpu sz) {
     struct cell_node *grid_cell = the_grid->first_cell;
 
-    real_cpu sx, sy, sz;
-    sx = 20000;
-    sy = 7000;
-    sz = 3000;
+    //real_cpu sx, sy, sz;
+    //sx = 20000;
+    //sy = 7000;
+    //sz = 3000;
 
     while(grid_cell != 0) {
         grid_cell->active = (grid_cell->center.x < sx) && (grid_cell->center.y < sy) && (grid_cell->center.z < sz);
@@ -1010,12 +1010,10 @@ int calc_num_refs(real_cpu start_h, real_cpu desired_h) {
     return num_refs;
 }
 
-void set_plain_fibrosis_source_sink_region (struct grid *the_grid, real_cpu phi, unsigned fib_seed, const double min_x, const double max_x, const double min_y,
-                                      const double max_y, const double min_z, const double max_z,
-                                      real_cpu source_sink_min_x, real_cpu source_sink_max_x, real_cpu side_length) {
-    log_info("Making %.2lf %% of cells inside the region inactive\n", phi * 100.0);
+void set_cuboid_sphere_fibrosis_with_conic_path(struct grid *the_grid, real_cpu phi, real_cpu plain_center, real_cpu sphere_radius, real_cpu bz_size, real_cpu bz_radius,
+                               unsigned fib_seed, real_cpu cone_slope) {
 
-    struct cell_node *grid_cell;
+    log_info("Making %.2lf %% of cells inactive\n", phi * 100.0f);
 
     if(fib_seed == 0)
         fib_seed = (unsigned)time(NULL) + getpid();
@@ -1024,30 +1022,63 @@ void set_plain_fibrosis_source_sink_region (struct grid *the_grid, real_cpu phi,
 
     log_info("Using %u as seed\n", fib_seed);
 
-    real_cpu a1 = (2.0*side_length) / (side_length - 2*source_sink_min_x);
-    real_cpu b1 = -source_sink_min_x*a1;
-    real_cpu a2 = (2.0*side_length) / (side_length - 2*source_sink_max_x);
-    real_cpu b2 = -source_sink_max_x*a2;
+    real_cpu bz_radius_2 = pow(bz_radius, 2.0);
+    real_cpu sphere_radius_2 = pow(sphere_radius, 2.0);
+    struct cell_node *grid_cell;
 
     grid_cell = the_grid->first_cell;
     while(grid_cell != 0) {
-        real center_x = grid_cell->center.x;
-        real center_y = grid_cell->center.y;
-        real center_z = grid_cell->center.z;
+        //Calcula distância da célula para o centro da malha
+        real_cpu distance = pow(grid_cell->center.x - plain_center, 2.0) + pow(grid_cell->center.y - plain_center, 2.0);
+        real_cpu h_distance = abs(grid_cell->center.x - plain_center); 
+        
+        if(grid_cell->active) {
 
-        if(center_x >= min_x && center_x <= max_x && center_y >= min_y && center_y <= max_y && center_z >= min_z && center_z <= max_z
-            && (center_y > a1*center_x + b1 || center_y > a2*center_x + b2)) {
-            if(grid_cell->active) {
-                real_cpu p = (real_cpu)(rand()) / (RAND_MAX);
-                if(p < phi) {
-                    grid_cell->active = false;
-                }
+            INITIALIZE_FIBROTIC_INFO(grid_cell);
 
-                INITIALIZE_FIBROTIC_INFO(grid_cell);
-                FIBROTIC(grid_cell) = true;
+            if(distance <= bz_radius_2) {
+                //Dentro da border zone 
+                if(distance <= sphere_radius_2) {
+                    //dentro da "esfera"
+                    if(h_distance < cone_slope*grid_cell->center.y*plain_center) //abs(cone_slope*grid_cell->center.y)
+                    {
+                        FIBROTIC(grid_cell) = true;
+                        
+                        //Dentro do cone
+                    }
+                    else{
+                        grid_cell->active = false;
+                        grid_cell->can_change = false;
+                        FIBROTIC(grid_cell) = true;
+                    }
+                } else {
+                    BORDER_ZONE(grid_cell) = true;
+                }
             }
         }
+        grid_cell = grid_cell->next;
+    }
+
+    grid_cell = the_grid->first_cell;
 
+    while(grid_cell != 0) {
+        if(grid_cell->active) {
+            if(FIBROTIC(grid_cell)) {
+                real_cpu p = (real_cpu)(rand()) / (RAND_MAX);
+                if(p < phi)
+                    grid_cell->active = false;
+                grid_cell->can_change = false;
+            } else if(BORDER_ZONE(grid_cell)) {
+                real_cpu distance_from_center = sqrt((grid_cell->center.x - plain_center) * (grid_cell->center.x - plain_center) +
+                                                     (grid_cell->center.y - plain_center) * (grid_cell->center.y - plain_center));
+                distance_from_center = (distance_from_center - sphere_radius) / bz_size;
+                real_cpu phi_local = phi - phi * distance_from_center;
+                real_cpu p = (real_cpu)(rand()) / (RAND_MAX);
+                if(p < phi_local)
+                    grid_cell->active = false;
+                grid_cell->can_change = false;
+            }
+        }
         grid_cell = grid_cell->next;
     }
-}
+}
diff --git a/src/domains_library/domain_helpers.h b/src/domains_library/domain_helpers.h
@@ -17,7 +17,7 @@ int set_cuboid_domain_mesh(struct grid *the_grid, real_cpu start_dx, real_cpu st
                            real_cpu side_length_z);
 int set_square_mesh(struct config *config, struct grid *the_grid);
 
-void set_benchmark_domain(struct grid *the_grid);
+void set_benchmark_domain(struct grid *the_grid, real_cpu sx, real_cpu sy, real_cpu sz);
 void set_cuboid_domain(struct grid *the_grid, real_cpu sizeX, real_cpu sizeY, real_cpu sizeZ);
 
 void set_custom_mesh(struct grid *the_grid, const char *file_name, size_t size, char *read_format);
@@ -56,6 +56,9 @@ uint32_t set_custom_mesh_from_file(struct grid *the_grid, const char *mesh_file,
 
 void set_cube_sphere_fibrosis(struct grid *the_grid, real_cpu phi, real_cpu sphere_center[3], real_cpu sphere_radius, unsigned fib_seed);
 
+void set_cuboid_sphere_fibrosis_with_conic_path (struct grid *the_grid, real_cpu phi, real_cpu plain_center, real_cpu sphere_radius, real_cpu bz_size, real_cpu bz_radius,
+                               unsigned fib_seed, real_cpu cone_slope);
+
 int calc_num_refs(real_cpu start_h, real_cpu desired_h);
 
 #endif // MONOALG3D_DOMAIN_HELPERS_H
