CGAL Polyline Simplification results in self-intersection

Question

I'm currently having a bit of trouble with CGAL's Polyline Simplification.

More specifically, for the following example, PS::simplify(ct, Cost(), Stop(0.2)) results in a self-intersecting polyline. In the image below, the blue polyline is the input polyline into PS::simplify() while the green polyline is the resulting (output) polyline. The red arrow points to the self-intersection in the resulting polyline.

Further below, I have copied and pasted my code from 2 files simplify_test.cpp and CMakeLists.txt. With the required libraries installed, to run this example, you may place them in the same directory, cd to that directory, and run the following in your terminal:

$ cmake .
$ make
$ ./simplify_test

This outputs the resulting coordinates (green polyline in the image below), which contain the self-intersection.

I would like to know:
(1) why this results in a self-intersection
and (2) what can be done so that it does not result in a self-intersection.

Thank you for your help!

Here is the simplification code in a file named simplify_test.cpp:

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Polyline_simplification_2/simplify.h>

#include <iostream>
#include <string>
#include <vector>

typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Point_2<K> Point;
typedef std::vector<Point> Polyline;
namespace PS = CGAL::Polyline_simplification_2;
typedef PS::Vertex_base_2<K> Vb;
typedef CGAL::Constrained_triangulation_face_base_2<K> Fb;
typedef CGAL::Triangulation_data_structure_2<Vb, Fb> TDS;
typedef CGAL::Constrained_Delaunay_triangulation_2<K, TDS, CGAL::Exact_predicates_tag> CDT;
typedef CGAL::Constrained_triangulation_plus_2<CDT> CT;
typedef PS::Stop_below_count_ratio_threshold Stop;
typedef PS::Squared_distance_cost Cost;

void print_coords(Polyline polyline) {
  std::cout << std::endl;
  std::cout << "Simplified coordinates:" << std::endl << std::endl;

  // Print out the simplified coordinates
  unsigned int i = 0;
  for (Point coord : polyline) {
    std::cout << "[ ";
    std::cout << coord.x();
    std::cout << ", ";
    std::cout << coord.y();
    std::cout << " ]";
    if (i != polyline.size() - 1) std::cout << ",";
    i++;
  }
  std::cout << std::endl << std::endl;
}

void simplify_test() {
  // Hard code a minimum working example where running PS::simplify results in
  // self-intersections. There are no self-intersections when {27, 9} is
  // omitted.
  std::vector<std::vector<int>> coords = {
      {64, 20}, {33, 27}, {27, 9}, {33, 18}, {44, 18}, {44, 8},
      {24, 0},  {0, 13},  {9, 49}, {84, 41}, {83, 29}, {64, 20},
  };
  // Simplification outputs:
  // [ 64, 20 ],[ 27, 9 ],[ 44, 18 ],[ 24, 0 ],
  // [ 9, 49 ],[ 83, 29 ],[ 64, 20 ],[ 64, 20 ]

  // Create polyline for simplifying later
  Polyline polyline;

  // Insert coordinates into polyline
  for (std::vector<int> coord : coords) {
    Point pt(coord[0], coord[1]);
    polyline.push_back(pt);
  }

  // Insert polyline into ct and run simplify()
  CT ct;
  ct.insert_constraint(polyline.begin(), polyline.end());
  PS::simplify(ct, Cost(), Stop(0.2));
  Polyline polyline_simplified;

  // Transfer simplified coordinates from ct to polyline for easy handling
  auto cit = ct.constraints_begin();
  for (auto vit = ct.points_in_constraint_begin(*cit);
       vit != ct.points_in_constraint_end(*cit); vit++) {
    polyline_simplified.push_back(*vit);
  }

  print_coords(polyline_simplified);
}

int main() {
  simplify_test();
}

Here is the CMakeLists.txt file.

cmake_minimum_required(VERSION 3.1)
project(simplify_test LANGUAGES CXX)

set(CMAKE_CXX_STANDARD 20)
set(CMAKE_BUILD_TYPE Release)

# Detecting appropriate compiler
if (APPLE)
  set(CMAKE_C_COMPILER "/usr/local/opt/llvm/bin/clang")
  set(CMAKE_CXX_COMPILER "/usr/local/opt/llvm/bin/clang++")
elseif(UNIX) # implicit AND NOT APPLE
  set(CMAKE_CXX_COMPILER "g++-10")
  set(CMAKE_C_COMPILER "gcc-10")
endif()

# Finding appropriate packages
find_package(CGAL)

# Adding executables needed
add_executable(simplify_test ./simplify_test.cpp)

# Linking appropriate libraries required
target_link_libraries(
  simplify_test
  CGAL::CGAL
)

Answer 1

I further reduced it to a polyline std::vector<std::vector<int> > coords = { {64, 20}, {33, 27}, {27, 9}, {33, 18}, {44, 18}, {24, 0} };

You found a bug in the points in constraint iterator. As a workaround you should use the Vertex_in_constraint_iterator and then call the method point(). I will file an issue on github and a fix.