diff --git a/Polygon/include/CGAL/Multipolygon_with_holes_2.h b/Polygon/include/CGAL/Multipolygon_with_holes_2.h
index bdfd3cb47bc1..9228850edf54 100644
--- a/Polygon/include/CGAL/Multipolygon_with_holes_2.h
+++ b/Polygon/include/CGAL/Multipolygon_with_holes_2.h
@@ -24,14 +24,14 @@ namespace CGAL {
 /*! \ingroup PkgPolygon2Ref
  *
  * The class `Multipolygon_with_holes_2` models the concept `MultipolygonWithHoles_2`.
- * It is parameterized with two types (`Kernel` and `Container`) that are used to instantiate
- * the types `Polygon_2<Kernel,Container>` and `Polygon_with_holes_2<Kernel,Container>`.
+ * It is parameterized with two types (`Kernel` and `Container_`) that are used to instantiate
+ * the types `Polygon_2<Kernel,Container_>` and `Polygon_with_holes_2<Kernel,Container_>`.
  * The latter is used to represent each polygon with holes. The former is converted to the latter.
  *
  * \cgalModels{MultipolygonWithHoles_2}
  */
 template <class Kernel,
-          class Container = std::vector<typename Kernel::Point_2>>
+          class Container_ = std::vector<typename Kernel::Point_2>>
 class Multipolygon_with_holes_2 {
 public:
   /// \name Definition
@@ -39,14 +39,15 @@ class Multipolygon_with_holes_2 {
   /// @{
 
   /// polygon type
-  using Polygon_2 = CGAL::Polygon_2<Kernel, Container>;
+  using Polygon_2 = CGAL::Polygon_2<Kernel, Container_>;
 
   /// polygon with holes type
-  using Polygon_with_holes_2 = CGAL::Polygon_with_holes_2<Kernel, Container>;
+  using Polygon_with_holes_2 = CGAL::Polygon_with_holes_2<Kernel, Container_>;
 
   /// @}
 
   using Traits = Kernel;
+  using Container = Container_;
   using value_type = Polygon_with_holes_2;
   using Polygon_with_holes_container = std::deque<Polygon_with_holes_2>;
 
@@ -183,10 +184,10 @@ order.
 
 \relates Multipolygon_with_holes_2
 */
-template <class Kernel, class Container>
+template <class Kernel, class Container_>
 std::ostream& operator<<(std::ostream& os,
-                         const Multipolygon_with_holes_2<Kernel, Container>& mp) {
-  typename Multipolygon_with_holes_2<Kernel, Container>::Polygon_with_holes_const_iterator i;
+                         const Multipolygon_with_holes_2<Kernel, Container_>& mp) {
+  typename Multipolygon_with_holes_2<Kernel, Container_>::Polygon_with_holes_const_iterator i;
 
   switch(IO::get_mode(os)) {
     case IO::ASCII :
@@ -214,11 +215,11 @@ std::ostream& operator<<(std::ostream& os,
   }
 }
 
-template <class Transformation, class Kernel, class Container>
-Multipolygon_with_holes_2<Kernel, Container> transform(const Transformation& t,
-                                                       const Multipolygon_with_holes_2<Kernel, Container>& mp)
+template <class Transformation, class Kernel, class Container_>
+Multipolygon_with_holes_2<Kernel, Container_> transform(const Transformation& t,
+                                                       const Multipolygon_with_holes_2<Kernel, Container_>& mp)
   {
-    Multipolygon_with_holes_2<Kernel, Container> result;
+    Multipolygon_with_holes_2<Kernel, Container_> result;
     for(const auto& pwh : mp.polygons_with_holes()){
       result.add_polygon_with_holes(transform(t, pwh));
     }
diff --git a/Polygon/include/CGAL/Polygon_2.h b/Polygon/include/CGAL/Polygon_2.h
index c3300d939fa4..cfe9b4f34723 100644
--- a/Polygon/include/CGAL/Polygon_2.h
+++ b/Polygon/include/CGAL/Polygon_2.h
@@ -60,8 +60,8 @@ namespace CGAL {
 /// algorithms were used and what complexity they have.
 ///
 
-template <class Traits_P, class Container_P
-        = std::vector<typename Traits_P::Point_2> >
+template <class Traits_, class Container_
+        = std::vector<typename Traits_::Point_2> >
 class Polygon_2 {
 
   public:
@@ -70,33 +70,33 @@ class Polygon_2 {
     /// @{
 
     /// The traits type.
-    typedef Traits_P Traits;
+    typedef Traits_ Traits;
     /// The container type.
-    typedef Container_P Container;
+    typedef Container_ Container;
 
     /// The number type of the coordinates of the points of the polygon.
-    typedef typename Traits_P::FT FT;
+    typedef typename Traits_::FT FT;
     /// The point type of the polygon.
-    typedef typename Traits_P::Point_2 Point_2;
+    typedef typename Traits_::Point_2 Point_2;
     /// The type of a segment between two points of the polygon.
-    typedef typename Traits_P::Segment_2 Segment_2;
+    typedef typename Traits_::Segment_2 Segment_2;
 
     /// @}
 
-    typedef typename Container_P::difference_type difference_type;
-    typedef typename Container_P::value_type value_type;
-    typedef typename Container_P::pointer pointer;
-    typedef typename Container_P::reference reference;
-    typedef typename Container_P::const_reference const_reference;
+    typedef typename Container_::difference_type difference_type;
+    typedef typename Container_::value_type value_type;
+    typedef typename Container_::pointer pointer;
+    typedef typename Container_::reference reference;
+    typedef typename Container_::const_reference const_reference;
 
 
     //-------------------------------------------------------//
     // this intermediary step is required by Sun C++ 4.1
-    typedef typename Container_P::iterator iterator;
-    typedef typename Container_P::const_iterator const_iterator;
+    typedef typename Container_::iterator iterator;
+    typedef typename Container_::const_iterator const_iterator;
     //-------------------------------------------------------//
     typedef typename Container::iterator       Vertex_const_iterator;
-    typedef Polygon_circulator<Container_P>    Vertex_const_circulator;
+    typedef Polygon_circulator<Container_>    Vertex_const_circulator;
 
     /// \name Iterators
     ///
@@ -141,11 +141,11 @@ class Polygon_2 {
   //
 #else
     typedef Vertex_const_circulator Vertex_circulator;
-    typedef Polygon_2_edge_iterator<Traits_P,Container_P> Edge_const_iterator;
-    typedef Polygon_2_const_edge_circulator<Traits_P,
-                                            Container_P> Edge_const_circulator;
+    typedef Polygon_2_edge_iterator<Traits_,Container_> Edge_const_iterator;
+    typedef Polygon_2_const_edge_circulator<Traits_,
+                                            Container_> Edge_const_circulator;
 
-    typedef Polygon_2_edge_iterator<Traits_P,Container_P,
+    typedef Polygon_2_edge_iterator<Traits_,Container_,
                                     Tag_false> Vertex_pair_iterator;
 
   typedef Iterator_range<Edge_const_iterator> Edges;
@@ -162,7 +162,7 @@ class Polygon_2 {
     Polygon_2(const Traits & p_traits) : traits(p_traits) {}
 
     // Move constructor
-    // Polygon_2(Polygon_2<Traits_P,Container_P>&& polygon) = default;
+    // Polygon_2(Polygon_2<Traits_,Container_>&& polygon) = default;
 
     /// Creates a polygon with vertices from the sequence
     /// defined by the range \c [first,last).
@@ -260,7 +260,7 @@ class Polygon_2 {
     {
       if (size() <= 1)
         return;
-      typename Container_P::iterator i = d_container.begin();
+      typename Container_::iterator i = d_container.begin();
       std::reverse(++i, d_container.end());
     }
 
@@ -500,32 +500,32 @@ class Polygon_2 {
       { return d_container.empty(); }
 
     /// Returns a const reference to the sequence of vertices of the polygon.
-    const Container_P& container() const
+    const Container_& container() const
       { return d_container; }
 
     /// Returns a reference to the sequence of vertices of the polygon.
-    Container_P& container()
+    Container_& container()
       { return d_container; }
 
     /// Returns an iterator to the first vertex of the polygon.
-    typename Container_P::iterator begin()
+    typename Container_::iterator begin()
     {
        return container().begin();
     }
     /// Returns an iterator to the element after the last vertex of the polygon.
-    typename Container_P::iterator end()
+    typename Container_::iterator end()
     {
        return container().end();
     }
 
     /// Returns a const iterator to the first vertex of the polygon.
-    const typename Container_P::const_iterator begin() const
+    const typename Container_::const_iterator begin() const
     {
        return container().begin();
     }
 
     /// Returns a const iterator to the element after the last vertex of the polygon.
-    const typename Container_P::const_iterator end() const
+    const typename Container_::const_iterator end() const
     {
        return container().end();
     }
@@ -544,14 +544,14 @@ class Polygon_2 {
 
     /// @}
 
-    bool identical(const Polygon_2<Traits_P,Container_P> &q) const
+    bool identical(const Polygon_2<Traits_,Container_> &q) const
       { return this == &q; }
 
-    Traits_P const &traits_member() const { return traits;}
+    Traits_ const &traits_member() const { return traits;}
 
   private:
-    Container_P d_container;
-    Traits_P traits;
+    Container_ d_container;
+    Traits_ traits;
 };
 
 
@@ -563,23 +563,23 @@ class Polygon_2 {
 /// equal to the vertices of `p1`. Note that the template argument
 /// `%Container` of `p1` and `p2` may be different.
 /// \memberof Polygon_2
-template <class Traits_P, class Container1_P, class Container2_P>
-bool operator==( const Polygon_2<Traits_P,Container1_P> &p1,
-                 const Polygon_2<Traits_P,Container2_P> &p2 );
+template <class Traits_, class Container1_P, class Container2_P>
+bool operator==( const Polygon_2<Traits_,Container1_P> &p1,
+                 const Polygon_2<Traits_,Container2_P> &p2 );
 
 /// Test for inequality.
 /// \memberof Polygon_2
-template <class Traits_P, class Container1_P, class Container2_P>
+template <class Traits_, class Container1_P, class Container2_P>
 inline
 bool
-operator!=(const Polygon_2<Traits_P,Container1_P> &p1,
-           const Polygon_2<Traits_P,Container2_P> &p2);
+operator!=(const Polygon_2<Traits_,Container1_P> &p1,
+           const Polygon_2<Traits_,Container2_P> &p2);
 
 /// Returns the image of the polygon \c p under the transformation \c t.
 /// \relates Polygon_2
-template <class Transformation, class Traits_P, class Container_P>
-Polygon_2<Traits_P,Container_P>
-transform(const Transformation& t, const Polygon_2<Traits_P,Container_P>& p);
+template <class Transformation, class Traits_, class Container_>
+Polygon_2<Traits_,Container_>
+transform(const Transformation& t, const Polygon_2<Traits_,Container_>& p);
 
 /// @} // global operators
 
@@ -591,14 +591,14 @@ transform(const Transformation& t, const Polygon_2<Traits_P,Container_P>& p);
 /// Reads a polygon from stream `is` and assigns it to `p`.
 /// \pre The extract operator must be defined for `Point_2`.
 /// \relates Polygon_2
-template <class Traits_P, class Container_P>
-std::istream &operator>>(std::istream &is, Polygon_2<Traits_P,Container_P>& p);
+template <class Traits_, class Container_>
+std::istream &operator>>(std::istream &is, Polygon_2<Traits_,Container_>& p);
 
 /// Inserts the polygon `p` into the stream `os`.
 /// \pre The insert operator must be defined for `Point_2`.
 /// \relates Polygon_2
-template <class Traits_P, class Container_P>
-std::ostream &operator<<(std::ostream &os, const Polygon_2<Traits_P,Container_P>& p);
+template <class Traits_, class Container_>
+std::ostream &operator<<(std::ostream &os, const Polygon_2<Traits_,Container_>& p);
 
 /// @} // IO
 
@@ -612,11 +612,11 @@ std::ostream &operator<<(std::ostream &os, const Polygon_2<Traits_P,Container_P>
 
 namespace CGAL {
 
-template <class Traits_P, class Container1_P, class Container2_P>
+template <class Traits_, class Container1_P, class Container2_P>
 inline
 bool
-operator!=(const Polygon_2<Traits_P,Container1_P> &x,
-           const Polygon_2<Traits_P,Container2_P> &y)
+operator!=(const Polygon_2<Traits_,Container1_P> &x,
+           const Polygon_2<Traits_,Container2_P> &y)
 {
   return !(x==y);
 }
diff --git a/Polygon/include/CGAL/Polygon_with_holes_2.h b/Polygon/include/CGAL/Polygon_with_holes_2.h
index e7178d650cae..7d581a281502 100644
--- a/Polygon/include/CGAL/Polygon_with_holes_2.h
+++ b/Polygon/include/CGAL/Polygon_with_holes_2.h
@@ -29,21 +29,22 @@ namespace CGAL {
 
 The class `Polygon_with_holes_2` models the concept `GeneralPolygonWithHoles_2`.
 It represents a linear polygon with holes. It is parameterized with two
-types (`Kernel` and `Container`) that are used to instantiate
-the type `Polygon_2<Kernel,Container>`. This polygon type is used to
+types (`Kernel` and `Container_`) that are used to instantiate
+the type `Polygon_2<Kernel,Container_>`. This polygon type is used to
 represent the outer boundary and the boundary of the holes (if any exist).
 
 \cgalModels{GeneralPolygonWithHoles_2}
 
 */
 template <class Kernel,
-          class Containter = std::vector<typename Kernel::Point_2> >
+          class Container_ = std::vector<typename Kernel::Point_2> >
 class Polygon_with_holes_2 :
-  public General_polygon_with_holes_2<CGAL::Polygon_2<Kernel, Containter> >
+  public General_polygon_with_holes_2<CGAL::Polygon_2<Kernel, Container_> >
 {
 public:
   typedef Kernel                                     Traits;
-  typedef CGAL::Polygon_2<Kernel, Containter>        Polygon_2;
+  typedef Container_                                 Container;
+  typedef CGAL::Polygon_2<Kernel, Container>         Polygon_2;
   typedef General_polygon_with_holes_2<Polygon_2>    Base;
   typedef typename Base::Hole_const_iterator         Hole_const_iterator;
   typedef typename Base::Size                        Size;
@@ -90,11 +91,11 @@ class Polygon_with_holes_2 :
 
 };
 
-  template <class Transformation, class Kernel, class Container>
-  Polygon_with_holes_2<Kernel,Container> transform(const Transformation& t,
-                                                   const Polygon_with_holes_2<Kernel,Container>& pwh)
+  template <class Transformation, class Kernel, class Container_>
+  Polygon_with_holes_2<Kernel,Container_> transform(const Transformation& t,
+                                                   const Polygon_with_holes_2<Kernel,Container_>& pwh)
   {
-    Polygon_with_holes_2<Kernel,Container> result(transform(t, pwh.outer_boundary()));
+    Polygon_with_holes_2<Kernel,Container_> result(transform(t, pwh.outer_boundary()));
     for(const auto& hole : pwh.holes()){
       result.add_hole(transform(t, hole));
     }
diff --git a/Polygon_repair/doc/Polygon_repair/PackageDescription.txt b/Polygon_repair/doc/Polygon_repair/PackageDescription.txt
index 31e905f74c2a..bd014eb61fd0 100644
--- a/Polygon_repair/doc/Polygon_repair/PackageDescription.txt
+++ b/Polygon_repair/doc/Polygon_repair/PackageDescription.txt
@@ -33,6 +33,8 @@ The %union and the %intersection rule enable to combine similar polygons. }
 
 \cgalCRPSection{Functions}
 - `CGAL::Polygon_repair::repair()`
+- `CGAL::Polygon_repair::join()`
+- `CGAL::Polygon_repair::intersect()`
 
 \cgalCRPSection{Repair Rules}
 - `CGAL::Polygon_repair::Even_odd_rule`
diff --git a/Polygon_repair/doc/Polygon_repair/Polygon_repair.txt b/Polygon_repair/doc/Polygon_repair/Polygon_repair.txt
index ec41c0d9ce56..e5dc97ba11b1 100644
--- a/Polygon_repair/doc/Polygon_repair/Polygon_repair.txt
+++ b/Polygon_repair/doc/Polygon_repair/Polygon_repair.txt
@@ -103,7 +103,7 @@ Tbd.
 
 \section SectionPolygonRepair_UnionIntersection Union and Intersection Rule
 
-Tbd.
+Given several valid polygons this rule computes their union or intersection.
 
 \section SubsectionPolygonRepair_Notes Notes on the Output
 
diff --git a/Polygon_repair/include/CGAL/Polygon_repair/Boolean.h b/Polygon_repair/include/CGAL/Polygon_repair/Boolean.h
index e2e76113741d..4ea0eb62f397 100644
--- a/Polygon_repair/include/CGAL/Polygon_repair/Boolean.h
+++ b/Polygon_repair/include/CGAL/Polygon_repair/Boolean.h
@@ -37,7 +37,7 @@ namespace Polygon_repair {
 \tparam Kernel must be
 */
 
-template <typename Kernel>
+template <typename Kernel, typename Container_>
 class Boolean {
 
 private:
@@ -60,10 +60,11 @@ class Boolean {
   };
 
   using K = Kernel;
+  using Container = Container_;
   using Point_2 = typename K::Point_2;
-  using Polygon_2 = CGAL::Polygon_2<K>;
-  using Polygon_with_holes_2 = CGAL::Polygon_with_holes_2<K>;
-  using Multipolygon_with_holes_2 = CGAL::Multipolygon_with_holes_2<K>;
+  using Polygon_2 = CGAL::Polygon_2<K,Container>;
+  using Polygon_with_holes_2 = CGAL::Polygon_with_holes_2<K,Container>;
+  using Multipolygon_with_holes_2 = CGAL::Multipolygon_with_holes_2<K,Container>;
 
   using Itag = std::conditional_t<std::is_floating_point_v<typename K::FT>, Exact_predicates_tag, Exact_intersections_tag>;
   using Vb = CGAL::Triangulation_vertex_base_2<K>;
diff --git a/Polygon_repair/include/CGAL/Polygon_repair/repair.h b/Polygon_repair/include/CGAL/Polygon_repair/repair.h
index a4b4a8745f13..f33368aa2e0b 100644
--- a/Polygon_repair/include/CGAL/Polygon_repair/repair.h
+++ b/Polygon_repair/include/CGAL/Polygon_repair/repair.h
@@ -44,9 +44,9 @@ class Polygon_repair;
 
 /// \ingroup PkgPolygonRepairFunctions
 /// repairs polygon `p` using the given rule
-/// \tparam Kernel parameter of the input and output polygons
+/// \tparam Kernel parameter of the input and output polygons. Must be model of `ConstrainedDelaunayTriangulationTraits_2 `
 /// \tparam Container parameter of the input and output polygons
-///  \tparam Rule must be `Even_odd_rule`
+///  \tparam Rule must be `Even_odd_rule` or `Non_zero_rule`
 template <class Kernel, class Container, class Rule = Even_odd_rule>
 Multipolygon_with_holes_2<Kernel, Container> repair(const Polygon_2<Kernel, Container>& p , Rule = Rule())
 {
@@ -61,9 +61,9 @@ Multipolygon_with_holes_2<Kernel, Container> repair(const Polygon_2<Kernel, Cont
 
 /// \ingroup PkgPolygonRepairFunctions
 /// repairs polygon with holes `p` using the given rule
-/// \tparam Kernel parameter of the input and output polygons
+/// \tparam Kernel parameter of the input and output polygons. Must be model of `ConstrainedDelaunayTriangulationTraits_2 `
 /// \tparam Container parameter of the input and output polygons
-///  \tparam Rule must be `Even_odd_rule`
+///  \tparam Rule must be `Even_odd_rule` or `Non_zero_rule`
 template <class Kernel, class Container, class Rule = Even_odd_rule>
 Multipolygon_with_holes_2<Kernel, Container> repair(const Polygon_with_holes_2<Kernel, Container>& p, Rule = Rule())
 {
@@ -90,9 +90,10 @@ Multipolygon_with_holes_2<Kernel, Container> repair(const Polygon_with_holes_2<K
 
 /// \ingroup PkgPolygonRepairFunctions
 /// repairs multipolygon with holes `p` using the given rule
-/// \tparam Kernel parameter of the input and output polygons
+/// \tparam Kernel parameter of the input and output polygons. Must be model of `ConstrainedDelaunayTriangulationTraits_2 `
 /// \tparam Container parameter of the input and output polygons
-///  \tparam Rule must be `Even_odd_rule` or `Non_zero_rule`
+///  \tparam Rule may be any rule
+/// \pre If the rule is the `Union_rule` or `Non_zero_rule`, each polygon with hole must be free of self-intersections
 template <class Kernel, class Container, class Rule = Even_odd_rule>
 Multipolygon_with_holes_2<Kernel, Container> repair(const Multipolygon_with_holes_2<Kernel, Container>& p, Rule = Rule())
 {
@@ -116,7 +117,7 @@ Multipolygon_with_holes_2<Kernel, Container> repair(const Multipolygon_with_hole
     }
   };
 
-  CGAL::Polygon_repair::Boolean<Kernel> bops;
+  CGAL::Polygon_repair::Boolean<Kernel,Container> bops;
   bops.insert(p);
   bops.mark_domains();
   Larger_than_zero ltz;
@@ -139,7 +140,7 @@ Multipolygon_with_holes_2<Kernel, Container> repair(const Multipolygon_with_hole
     }
   };
 
- CGAL::Polygon_repair::Boolean<Kernel> bops;
+ CGAL::Polygon_repair::Boolean<Kernel,Container> bops;
   bops.insert(p);
   bops.mark_domains();
   Equal equal(p.number_of_polygons_with_holes());
@@ -787,12 +788,13 @@ class Polygon_repair {
 
 
 /// \ingroup PkgPolygonRepairFunctions
-/// omputes the union of all polygons with holes in `p`
-/// \tparam K parameter of the input and output polygons
+/// computes the union of all polygons with holes in `p`
+/// \tparam Kernel parameter of the input and output polygons. Must be model of `ConstrainedDelaunayTriangulationTraits_2 `
+/// \tparam Container parameter of the input and output polygons
 /// \pre Each polygon with hole must be free of self-intersections
-template <typename K>
-Multipolygon_with_holes_2<K>
-join(const Multipolygon_with_holes_2<K>& pA)
+template <typename Kernel, typename Container>
+Multipolygon_with_holes_2<Kernel,Container>
+join(const Multipolygon_with_holes_2<Kernel,Container>& pa)
 {
   struct Larger_than_zero {
     bool operator()(int i) const
@@ -801,8 +803,8 @@ join(const Multipolygon_with_holes_2<K>& pA)
     }
   };
 
-  CGAL::Polygon_repair::Boolean<K> bops;
-  bops.insert(pA);
+  CGAL::Polygon_repair::Boolean<Kernel> bops;
+  bops.insert(pa);
   bops.mark_domains();
   Larger_than_zero ltz;
   return bops(ltz);
@@ -810,15 +812,21 @@ join(const Multipolygon_with_holes_2<K>& pA)
 
 /// \ingroup PkgPolygonRepairFunctions
 /// computes the union of two polygons
-/// \tparam K parameter of the output polygons
-/// \tparam PA must be `Polygon_2<K>`, or `Polygon_with_holes_2<K>`, or `Multipolygon_with_holes_2<K>`
-/// \tparam PB must be `Polygon_2<K>`, or `Polygon_with_holes_2<K>`, or `Multipolygon_with_holes_2<K>`
+/// \tparam Kernel parameter of the output polygons. Must be model of `ConstrainedDelaunayTriangulationTraits_2 `
+/// \tparam Container parameter of the input and output polygons
+/// \tparam PA must be `Polygon_2<Kernel, Container>`, or `Polygon_with_holes_2<Kernel, Container>`, or `Multipolygon_with_holes_2<Kernel, Container>`
+/// \tparam PB must be `Polygon_2<Kernel, Container>`, or `Polygon_with_holes_2<Kernel, Container>`, or `Multipolygon_with_holes_2<Kernel, Container>`
 /// \pre The polygons `pA` and `pB` must be free of self-intersections
-template <typename PA, typename PB, typename K = Default>
-decltype(auto) // Multipolygon_with_holes_2<K>
-join(const PA& pa, const PB& pb, const K& = Default())
+template <typename PA, typename PB, typename Kernel = Default, typename Container = Default>
+#ifdef DOXYGEN_RUNNING
+Multipolygon_with_holes_2<Kernel,Container>
+#else
+decltype(auto)
+#endif
+join(const PA& pa, const PB& pb, const Kernel& = Default(), const Container& = Default())
 {
-  typedef typename Default::Get<K, typename PA::Traits>::type Traits;
+  typedef typename Default::Get<Kernel, typename PA::Traits>::type Traits;
+  typedef typename Default::Get<Container, typename PA::Container>::type Container;
 
   struct Larger_than_zero {
     bool operator()(int i) const
@@ -827,7 +835,7 @@ join(const PA& pa, const PB& pb, const K& = Default())
     }
   };
 
-  CGAL::Polygon_repair::Boolean<Traits> bops;
+  CGAL::Polygon_repair::Boolean<Traits,Container> bops;
   bops.insert(pa);
   bops.insert(pb);
   bops.mark_domains();
@@ -838,11 +846,12 @@ join(const PA& pa, const PB& pb, const K& = Default())
 
 /// \ingroup PkgPolygonRepairFunctions
 /// computes the intersection of all polygons with holes in `p`
-/// \tparam K parameter of the input and output polygons
+/// \tparam Kernel parameter of the input and output polygons. Must be model of `ConstrainedDelaunayTriangulationTraits_2 `
+/// \tparam Container parameter of the input and output polygons
 /// \pre Each polygon with hole must be free of self-intersections
-template <typename K>
-Multipolygon_with_holes_2<K>
-intersect(const Multipolygon_with_holes_2<K>& p)
+template <typename Kernel, typename Container>
+Multipolygon_with_holes_2<Kernel,Container>
+intersect(const Multipolygon_with_holes_2<Kernel,Container>& p)
 {
   struct Equal  {
     int val;
@@ -856,7 +865,7 @@ intersect(const Multipolygon_with_holes_2<K>& p)
     }
   };
 
-  CGAL::Polygon_repair::Boolean<K> bops;
+  CGAL::Polygon_repair::Boolean<Kernel,Container> bops;
   bops.insert(p);
   bops.mark_domains();
   Equal equal(p.number_of_polygons_with_holes());
@@ -866,15 +875,21 @@ intersect(const Multipolygon_with_holes_2<K>& p)
 
 /// \ingroup PkgPolygonRepairFunctions
 /// Computes the intersection of two polygons
-/// \tparam K parameter of the output polygon
-/// \tparam PA must be `Polygon_2<K>`, or `Polygon_with_holes_2<K>`, or `Multipolygon_with_holes_2<K>`
-/// \tparam PB must be `Polygon_2<K>`, or `Polygon_with_holes_2<K>`, or `Multipolygon_with_holes_2<K>`
+/// \tparam Kernel parameter of the output polygon. Must be model of `ConstrainedDelaunayTriangulationTraits_2 `
+/// \tparam Container parameter of the input and output polygons
+/// \tparam PA must be `Polygon_2<Kernel, Container>`, or `Polygon_with_holes_2<Kernel, Container>`, or `Multipolygon_with_holes_2<Kernel, Container>`
+/// \tparam PB must be `Polygon_2<Kernel, Container>`, or `Polygon_with_holes_2<Kernel, Container>`, or `Multipolygon_with_holes_2<Kernel, Container>`
 /// \pre The polygons `pA` and `pB` must be free of self-intersections
-template <typename PA, typename PB, typename K = Default>
-decltype(auto)  // Multipolygon_with_holes_2<K>
-intersect(const PA& pa, const PB& pb, const K& = Default())
+template <typename PA, typename PB, typename Kernel = Default, typename Container = Default>
+#ifdef DOXYGEN_RUNNING
+Multipolygon_with_holes_2<Kernel, Container>
+#else
+decltype(auto)
+#endif
+intersect(const PA& pa, const PB& pb, const Kernel& = Default(), const Container& = Default())
 {
-  typedef typename Default::Get<K, typename PA::Traits>::type Traits;
+  typedef typename Default::Get<Kernel, typename PA::Traits>::type Traits;
+  typedef typename Default::Get<Container, typename PA::Container>::type Container;
 
   struct Equal  {
     bool operator()(int i) const
@@ -883,7 +898,7 @@ intersect(const PA& pa, const PB& pb, const K& = Default())
     }
   };
 
-  CGAL::Polygon_repair::Boolean<Traits> bops;
+  CGAL::Polygon_repair::Boolean<Traits,Container> bops;
   bops.insert(pa);
   bops.insert(pb);
   bops.mark_domains();