Add --timeout & --do_not_fail_on_intersections flag

.github/scripts/benchmark_maps.sh

@@ -9,7 +9,6 @@ BASE_BIN=$2
 OUT_JSON=$3
 
 MAP_ROOT="sample_data"
-WARMUP=1
 MIN_RUNS=6
 MAX_RUNS=15
 
@@ -81,10 +80,10 @@ for dir in "$MAP_ROOT"/*; do
     echo "Benchmarking:  $map_name"
 
     hyperfine --ignore-failure \
-      --warmup "$WARMUP" \
       --min-runs "$MIN_RUNS" \
       --max-runs "$MAX_RUNS" \
       --export-json "hf.json" \
+      --prepare 'bash -lc "if [ -e /proc/sys/vm/drop_caches ]; then sudo -n sh -c '\''sync; echo 3 > /proc/sys/vm/drop_caches'\'' >/dev/null 2>&1 || true; fi"' \
       --command-name main \
       "bash -c '[[ \$4 =~ world ]] && set -- \"\$1\" \"\$2\" \"\$3\" --world || set -- \"\$1\" \"\$2\" \"\$3\"; exec \"\$@\"' _ \
          \"$BASE_BIN\" \"$geo\" \"$csv\" \"$dir\"" \

include/inset_state.hpp

@@ -116,6 +116,8 @@ class InsetState
   // Whether convergence has been reached
   mutable bool converge_{true};
 
+  mutable bool intersections_found_{false};
+
   // Make default constructor private so that only
   // InsetState(const std::string, Arguments) can be called as constructor
   InsetState();

include/parse_arguments.hpp

@@ -35,7 +35,7 @@ struct Arguments {
   bool qtdt_method;
 
   // Should the polygons be simplified and densified?
-  bool simplify;
+  bool disable_simplification_densification;
 
   // If `rays` is true, we use the ray-shooting method to fill the grid cells.
   bool rays;
@@ -72,7 +72,11 @@ struct Arguments {
 
   bool verbose;
 
-  bool disable_triangulation_optimisation;
+  // Timeout in seconds
+  unsigned int timeout_in_seconds;
+
+  // Whether to exit gracefully if intersections are found
+  bool do_not_fail_on_intersections;
 
   // Column names in provided visual variables file (CSV)
   std::optional<std::string> id_col;

include/time_tracker.hpp

@@ -13,7 +13,12 @@ class TimeTracker
   std::unordered_map<std::string, std::chrono::milliseconds> durations_;
   std::string name_;
 
+  std::chrono::steady_clock::time_point program_start_;
+
 public:
+  TimeTracker();
+  explicit TimeTracker(std::string name);
+
   void set_name(std::string);
   void start(const std::string &task_name);
   void stop(const std::string &task_name);
@@ -22,6 +27,9 @@ class TimeTracker
 
   // Find the duration of a particular task
   std::chrono::milliseconds duration(const std::string &task_name) const;
+
+  // Total elapsed time from the CTOR of the object till now in seconds
+  double total_elapsed_time_in_seconds() const;
 };
 
 #endif  // TIME_TRACKER_H

run_maps.sh

@@ -46,7 +46,7 @@ for map_dir in */; do
       extra="--world"
     fi
 
-    command="../build/Release/cartogram \"../sample_data/${geojson_file}\" \"../sample_data/${csv_file}\" $extra"
+    command="../build/Debug/cartogram \"../sample_data/${geojson_file}\" \"../sample_data/${csv_file}\" $extra"
     echo "Command: $command"
 
     # Show the full path from the script location for the current map folder

src/cartogram_info/read_csv.cpp

@@ -221,7 +221,7 @@ static void check_validity_of_csv_ids(
 
   for (const auto &id : initial_id_order) {
     if (std::find(csv_ids.begin(), csv_ids.end(), id) == csv_ids.end()) {
-      std::cerr << "Warning: ID " << id << " in GeoJSON is not in CSV"
+      std::cerr << "WARNING: ID " << id << " in GeoJSON is not in CSV"
                 << std::endl;
       csv_data[id] =
         {{"area", "NA"}, {"color", ""}, {"label", ""}, {"inset_pos", "C"}};

src/inset_state/check_topology.cpp

@@ -30,13 +30,21 @@ void InsetState::holes_inside_polygons() const
 
 void InsetState::is_simple(const char *caller_func) const
 {
-  if (!args_.simplify)
+  if (args_.disable_simplification_densification)
     return;
 
   // Only check topology if simplification and densification is enabled.
   for (const auto &gd : geo_divs_) {
     for (const auto &pwh : gd.polygons_with_holes()) {
       if (!pwh.outer_boundary().is_simple()) {
+        intersections_found_ = true;
+        if (args_.do_not_fail_on_intersections) {
+          std::cerr << "WARNING: Outer boundary is not simple for GeoDiv "
+                    << gd.id();
+          std::cerr << ". is_simple() called from " << caller_func
+                    << std::endl;
+          continue;
+        }
         not_simple_polygon_ = pwh.outer_boundary();
         std::cerr << "ERROR: Outer boundary is not simple for GeoDiv "
                   << gd.id();
@@ -47,8 +55,16 @@ void InsetState::is_simple(const char *caller_func) const
           false);
         exit(1);
       }
+
       for (const auto &h : pwh.holes()) {
         if (!h.is_simple()) {
+          intersections_found_ = true;
+          if (args_.do_not_fail_on_intersections) {
+            std::cerr << "WARNING: Hole is not simple for GeoDiv " << gd.id();
+            std::cerr << ". is_simple() called from " << caller_func
+                      << std::endl;
+            continue;
+          }
           not_simple_polygon_ = h;
           std::cerr << "ERROR: Hole is not simple for GeoDiv " << gd.id();
           std::cerr << ". is_simple() called from " << caller_func

src/inset_state/inset_state.cpp

@@ -84,61 +84,6 @@ double InsetState::blur_width() const
   return blur_width;
 }
 
-bool InsetState::continue_integrating() const
-{
-
-  // Calculate all the necessary information to decide whether to continue
-  auto [max_area_err, worst_gd] = max_area_error();
-
-  // A GeoDiv is still above our area error threshold
-  bool area_error_above_threshold =
-    max_area_err > args_.max_permitted_area_error;
-
-  // Area expansion factor is above our threshold
-  // i.e. cartogram has become too big or too small
-  double area_drift = area_expansion_factor() - 1.0;
-  bool area_expansion_factor_above_threshold =
-    std::abs(area_drift) > max_permitted_area_drift;
-
-  // If both the above metrics are above our threshold
-  bool has_converged =
-    !area_error_above_threshold && !area_expansion_factor_above_threshold;
-
-  // Make sure to not continue endlesslely: cap at max_integrations
-  bool within_integration_limit = n_finished_integrations() < max_integrations;
-
-  // We continue if we are within the integration limit and have not converged
-  bool continue_integration =
-    (within_integration_limit && !has_converged) ||
-    (n_finished_integrations_ < args_.min_integrations);
-
-  // Actually hasn't converged, just reached integration limit
-  if (!within_integration_limit && !has_converged) {
-    converge_ = false;
-    std::cerr << "ERROR: Could not converge!" << std::endl;
-    if (area_error_above_threshold)
-      std::cerr << "Max area error above threshold!" << std::endl;
-    if (area_expansion_factor_above_threshold)
-      std::cerr << "Area expansion factor above threshold!" << std::endl;
-  }
-
-  // Print control output (at end of previous integration)
-  std::cerr << "Max. area err: " << max_area_err << ", GeoDiv: " << worst_gd
-            << std::endl;
-  std::cerr << "Current Area: " << geo_div_at_id(worst_gd).area()
-            << ", Target Area: " << target_area_at(worst_gd) << std::endl;
-  std::cerr << "Area drift: " << area_drift * 100.0 << "%" << std::endl;
-
-  if (continue_integration) {
-    // Print next integration information.
-    std::cerr << "\nIntegration number " << n_finished_integrations()
-              << std::endl;
-    std::cerr << "Dimensions : " << lx_ << " " << ly_ << std::endl;
-    std::cerr << "Number of Points: " << n_points() << std::endl;
-  }
-  return continue_integration;
-}
-
 Color InsetState::color_at(const std::string &id) const
 {
   try {

src/inset_state/integrate.cpp

@@ -17,7 +17,7 @@ void InsetState::preprocess()
     store_original_geo_divs();
   }
 
-  if (args_.simplify) {
+  if (!args_.disable_simplification_densification) {
     std::cerr << "Start of initial simplification of " << pos_ << std::endl;
 
     // Simplification reduces the number of points used to represent the
@@ -78,6 +78,73 @@ void InsetState::cleanup_after_integration()
   ref_to_fluxy_init().free();
 }
 
+bool InsetState::continue_integrating() const
+{
+  if (intersections_found_) {
+    std::cerr << "WARNING: Polygons do not remain simple during integration! "
+                 "Exiting gracefully."
+              << std::endl;
+    return false;
+  }
+
+  // Calculate all the necessary information to decide whether to continue
+  auto [max_area_err, worst_gd] = max_area_error();
+
+  // A GeoDiv is still above our area error threshold
+  bool area_error_above_threshold =
+    max_area_err > args_.max_permitted_area_error;
+
+  // Area expansion factor is above our threshold
+  // i.e. cartogram has become too big or too small
+  double area_drift = area_expansion_factor() - 1.0;
+  bool area_expansion_factor_above_threshold =
+    std::abs(area_drift) > max_permitted_area_drift;
+
+  // If both the above metrics are above our threshold
+  bool has_converged =
+    !area_error_above_threshold && !area_expansion_factor_above_threshold;
+
+  // Make sure to not continue endlesslely: cap at max_integrations
+  bool within_integration_limit = n_finished_integrations() < max_integrations;
+
+  // We continue if we are within the integration limit and have not converged
+  bool continue_integration =
+    (within_integration_limit && !has_converged) ||
+    (n_finished_integrations_ < args_.min_integrations);
+
+  // Actually hasn't converged, just reached integration limit
+  if (!within_integration_limit && !has_converged) {
+    converge_ = false;
+    std::cerr << "ERROR: Could not converge!" << std::endl;
+    if (area_error_above_threshold)
+      std::cerr << "Max area error above threshold!" << std::endl;
+    if (area_expansion_factor_above_threshold)
+      std::cerr << "Area expansion factor above threshold!" << std::endl;
+  }
+
+  // Print control output (at end of previous integration)
+  std::cerr << "Max. area err: " << max_area_err << ", GeoDiv: " << worst_gd
+            << std::endl;
+  std::cerr << "Current Area: " << geo_div_at_id(worst_gd).area()
+            << ", Target Area: " << target_area_at(worst_gd) << std::endl;
+  std::cerr << "Area drift: " << area_drift * 100.0 << "%" << std::endl;
+
+  if (timer.total_elapsed_time_in_seconds() > args_.timeout_in_seconds) {
+    std::cerr << "WARNING: Timeout of " << args_.timeout_in_seconds
+              << " seconds reached!" << std::endl;
+    return false;
+  }
+
+  if (continue_integration) {
+    // Print next integration information.
+    std::cerr << "\nIntegration number " << n_finished_integrations()
+              << std::endl;
+    std::cerr << "Dimensions : " << lx_ << " " << ly_ << std::endl;
+    std::cerr << "Number of Points: " << n_points() << std::endl;
+  }
+  return continue_integration;
+}
+
 void InsetState::integrate(ProgressTracker &progress_tracker)
 {
 

src/inset_state/project.cpp

@@ -8,7 +8,7 @@ bool InsetState::project()
   if (!create_delaunay_t())  // Triangle has flipped during triangulation
     return false;
 
-  if (args_.simplify) {
+  if (!args_.disable_simplification_densification) {
     densify_geo_divs_using_delaunay_t();
   }
 
@@ -28,7 +28,7 @@ bool InsetState::project()
       true);
   }
 
-  if (args_.simplify) {
+  if (!args_.disable_simplification_densification) {
 
     simplify(args_.target_points_per_inset);
   }

src/inset_state/simplify_inset.cpp

@@ -16,6 +16,13 @@ using Cost = CGAL::Polyline_simplification_2::Squared_distance_cost;
 
 void InsetState::simplify(const unsigned int target_points_per_inset)
 {
+  if (intersections_found_) {
+    std::cerr
+      << "WARNING: Intersections found previously; skipping simplification!"
+      << std::endl;
+    return;
+  }
+
   timer.start("Simplification");
   const size_t n_pts_before = n_points();
   std::cerr << n_pts_before << " points in inset. ";

src/misc/parse_arguments.cpp

@@ -1,13 +1,14 @@
 #include "parse_arguments.hpp"
 #include "constants.hpp"
+#include <limits>
 
 Arguments parse_arguments(const int argc, const char *argv[])
 {
   Arguments args;
 
   // Create parser for arguments using argparse.
   // From https://github.com/p-ranav/argparse
-  argparse::ArgumentParser arguments("./cartogram", "25.7");
+  argparse::ArgumentParser arguments("./cartogram", "25.8");
 
   // Positional argument accepting geometry file (GeoJSON, JSON) as input
   arguments.add_argument("geometry_file")
@@ -28,6 +29,12 @@ Arguments parse_arguments(const int argc, const char *argv[])
       "Integer: Number of starting grid cells along longer Cartesian "
       "coordinate axis");
 
+  arguments.add_argument("-T", "--timeout")
+    .default_value(std::numeric_limits<unsigned int>::max())
+    .scan<'u', unsigned int>()
+    .help(
+      "Integer: Maximum time (in seconds) allowed for cartogram generation");
+
   arguments.add_argument("-N", "--max_allowed_autoscale_grid_length")
     .default_value(max_allowed_autoscale_grid_length)
     .scan<'u', unsigned int>()
@@ -71,10 +78,6 @@ Arguments parse_arguments(const int argc, const char *argv[])
       "polygons")
     .default_value(false)
     .implicit_value(true);
-  arguments.add_argument("--disable_triangulation_optimisation")
-    .help("Boolean: Disable optimisation of maximum angle of triangulation")
-    .default_value(false)
-    .implicit_value(true);
   arguments.add_argument("--skip_projection")
     .help("Boolean: Skip projection to equal area")
     .default_value(false)
@@ -93,6 +96,14 @@ Arguments parse_arguments(const int argc, const char *argv[])
     // Currently, this help message is not accurate.
     .default_value(false)
     .implicit_value(true);
+
+  arguments.add_argument("--do_not_fail_on_intersections")
+    .help(
+      "Boolean: Whether to still produce cartogram if polygons do not remain "
+      "simple")
+    .default_value(false)
+    .implicit_value(true);
+
   arguments.add_argument("--output_shifted_insets")
     .help(
       "Boolean: Output repositioned insets in cartesian coordinates GeoJSON")
@@ -190,9 +201,8 @@ Arguments parse_arguments(const int argc, const char *argv[])
 
   // Set boolean values
   args.world = arguments.get<bool>("--world");
-  args.simplify = !arguments.get<bool>("--disable_simplify_and_densify");
-  args.disable_triangulation_optimisation =
-    arguments.get<bool>("--disable_triangulation_optimisation");
+  args.disable_simplification_densification =
+    arguments.get<bool>("--disable_simplify_and_densify");
   args.remove_tiny_polygons = arguments.get<bool>("--remove_tiny_polygons");
   args.min_polygon_area = arguments.get<double>("--minimum_polygon_area");
   args.max_permitted_area_error =
@@ -211,6 +221,10 @@ Arguments parse_arguments(const int argc, const char *argv[])
   args.plot_polygons = arguments.get<bool>("--plot_polygons");
   args.plot_quadtree = arguments.get<bool>("--plot_quadtree");
   args.output_shifted_insets = arguments.get<bool>("--output_shifted_insets");
+  args.do_not_fail_on_intersections =
+    arguments.get<bool>("--do_not_fail_on_intersections");
+
+  args.timeout_in_seconds = arguments.get<unsigned int>("--timeout");
 
   // arguments.present returns an optional
   args.id_col = arguments.present<std::string>("--id");
@@ -223,7 +237,7 @@ Arguments parse_arguments(const int argc, const char *argv[])
   args.verbose = arguments.get<bool>("--verbose");
 
   // Check whether n_points is specified but --simplify_and_densify not passed
-  if (!args.simplify) {
+  if (args.disable_simplification_densification) {
     std::cerr << "WARNING: Simplification and densification disabled! "
               << "Polygons will not simplified (or densified). "
               << "This may result and in polygon intersections. "