[pre-commit.ci] pre-commit autoupdate

.pre-commit-config.yaml

@@ -8,7 +8,7 @@ ci:
 
 repos:
   - repo: https://github.com/lorenzwalthert/precommit
-    rev: v0.4.3.9029
+    rev: v0.4.3.9030
     hooks:
       - id: parsable-R
       - id: no-browser-statement
@@ -49,19 +49,19 @@ repos:
         exclude: '(\.Rd|python/doc/source/reference/.*|test-doctest-.*)'
 
   - repo: https://github.com/tox-dev/tox-ini-fmt
-    rev: 1.7.2
+    rev: 1.9.0
     hooks:
       - id: tox-ini-fmt
 
   - repo: https://github.com/tox-dev/pyproject-fmt
-    rev: v2.25.2
+    rev: v2.29.4
     hooks:
       - id: pyproject-fmt
         additional_dependencies: ["tox>=4.12.1"]
 
   - repo: https://github.com/astral-sh/ruff-pre-commit
     # Ruff version.
-    rev: v0.15.21
+    rev: v0.16.7
     hooks:
       # Run the formatter.
       - id: ruff-format

c/NEWS.md

@@ -3,7 +3,11 @@
 ## 0.20
 
  * radixsort.h, insort.h: New.
-
+ * treap.h: New.
+ * pareto.c (pareto_rank_3d): Replace AVL-tree with Treap.
+ * nondominated.h (find_nondominated_3d_impl_sorted): Likewise.
+ * nondominated_kung.h (kung_merge_dim3): Likewise.
+ * hv3d_priv.h (hv3d_preprocessing): Likewise.
 
 ## 0.19.2
 

c/hv3d_priv.h

@@ -1,6 +1,9 @@
 #ifndef _HV3D_PRIV_H
 #define _HV3D_PRIV_H
 
+#define USE_AVL 0
+
+#if USE_AVL == 1
 typedef const double avl_item_t;
 typedef struct avl_node_t {
     struct avl_node_t *next;
@@ -29,6 +32,29 @@ new_avl_node(dlnode_t * restrict p, avl_node_t * restrict node)
     return node;
 }
 
+#else
+typedef struct TreapItem {
+    const double * z;
+    dlnode_t * dlnode;
+} TreapItem;
+
+static inline const double *
+treap_item_get_z(TreapItem item)
+{
+    return item.z;
+}
+
+#include "treap.h"
+
+static inline void
+hv3d_init_treap_node(TreapNode * restrict node, dlnode_t * restrict p)
+{
+    // Shift by -1 so that the split dimension is [0].
+    treap_node_init(node, p->x[1], (TreapItem) { .z = p->x - 1, .dlnode = p });
+}
+
+#endif
+
 /* Used by hvc3d.c and hv3dplus.c.
 
   This implements a variant of the 3D dimension-sweep algorithm by H. T. Kung,
@@ -61,6 +87,7 @@ hv3d_preprocessing(dlnode_t * restrict list, size_t n)
     assert(list+1 == list->next[0]);
     assert(list+2 == list->prev[0]);
 
+#if USE_AVL == 1
     avl_tree_t tree;
     avl_init_tree(&tree, qsort_cmp_pdouble_asc_y_des_x_nonzero);
     avl_node_t * tnodes = malloc((n+2) * sizeof(*tnodes));
@@ -131,6 +158,67 @@ hv3d_preprocessing(dlnode_t * restrict list, size_t n)
         }
         p = p->next[0];
     }
+#else
+    TreapNode *tnodes = malloc((n+2) * sizeof(*tnodes));
+    assert(tnodes != NULL);
+    // At the top we insert the first point, which is never dominated.
+    dlnode_t * p = (list+1)->next[0];
+    hv3d_init_treap_node(tnodes, p);
+    Treap tree;
+    treap_init_with_single_node(&tree, tnodes);
+    set_delimiters(p, list+1, list);
+
+    // After the top node, we insert sentinel 1 (-INF, ref[1])
+    hv3d_init_treap_node(tnodes + 1, list);
+    tnodes->right = tnodes + 1;
+    tnodes->right->priority = 0; // Push to the bottom
+    // Before the top node, we insert sentinel 2 (ref[0], -INF)
+    hv3d_init_treap_node(tnodes + 2, list + 1);
+    tnodes->left = tnodes + 2;
+    tnodes->left->priority = 0; // Push to the bottom
+    set_delimiters(p, tnodes->left->item.dlnode, tnodes->right->item.dlnode);
+
+    TreapNode * node = tnodes + 3;
+    const dlnode_t * stop = list+2;
+    _attr_maybe_unused dlnode_t * prev_p = p;
+    double pk0 = p->x[0], pk1 = p->x[1], _attr_maybe_unused pk2 = p->x[2];
+    for (p = p->next[0]; p != stop; p = p->next[0]) {
+        const double pj0 = p->x[0], pj1 = p->x[1], pj2 = p->x[2];
+        if ((pk0 > pj0) | (pk1 > pj1)) {
+            TreapNode *pred = treap_find_le(&tree, pj1);
+            assert(pred != NULL);
+            const double * prev_x = pred->item.dlnode->x;
+            if (prev_x[0] <= pj0) {
+                // pj is dominated by a point in the tree.
+#ifdef HVC_ONLY
+                if (all_equal_double(prev_x, p->x, 3))
+                    pred->item.dlnode->ignore = true; // It will have zero hvc.
+#endif
+                remove_from_z(p);
+                continue;
+            }
+            // pj is NOT dominated
+            hv3d_init_treap_node(node, p);
+            TreapNode * prev, *next;
+            (void) treap_insert_and_displace_get_bounds(&tree, node, &prev, &next);
+            node++;
+            // Check if the data structure is properly setup
+            assert(prev->item.dlnode->x[0] > pj0 && prev->item.dlnode->x[1] < pj1);
+            assert(next->item.dlnode->x[0] < pj0 && next->item.dlnode->x[1] > pj1);
+            set_delimiters(p, prev->item.dlnode, next->item.dlnode);
+
+            pk0 = pj0; pk1 = pj1; pk2 = pj2;
+            prev_p = p;
+        } else {
+            // pj is dominated by a previous point.
+#ifdef HVC_ONLY
+            if (pk0 == pj0 && pk1 == pj1 && pk2 == pj2)
+                prev_p->ignore = true; // It will have zero hvc.
+#endif
+            remove_from_z(p);
+        }
+    }
+#endif
     free(tnodes);
 #undef set_delimiters
 }

c/nondominated.h

@@ -1,12 +1,29 @@
 #ifndef NONDOMINATED_H
 #define NONDOMINATED_H
+/*****************************************************************************
+
+ Various algorithm for filtering dominated solutions
+
+ ---------------------------------------------------------------------
+
+ Copyright (C) 2026
+ Manuel Lopez-Ibanez <manuel.lopez-ibanez@manchester.ac.uk>
+
+ This Source Code Form is subject to the terms of the Mozilla Public
+ License, v. 2.0. If a copy of the MPL was not distributed with this
+ file, You can obtain one at https://mozilla.org/MPL/2.0/.
+
+*****************************************************************************/
 
 #include "config.h"
 #include <string.h> // memcpy
 #include <math.h> // INFINITY
 #include "sort.h"
 #include "radixsort.h"
 
+#define USE_AVL 0
+
+#if USE_AVL == 1
 typedef const double avl_item_t;
 typedef struct avl_node_t {
     struct avl_node_t *next;
@@ -19,6 +36,16 @@ typedef struct avl_node_t {
 } avl_node_t;
 
 #include "avl_tiny.h"
+#else
+typedef const double * TreapItem;
+static inline const double *
+treap_item_get_z(TreapItem item)
+{
+    return item;
+}
+
+#include "treap.h"
+#endif
 
 enum objs_agree_t { AGREE_MINIMISE = -1, AGREE_NONE = 0, AGREE_MAXIMISE = 1 };
 
@@ -335,6 +362,7 @@ find_nondominated_3d_impl_sorted(const double ** restrict rows, size_t size,
                                  const bool find_dominated)
 {
     ASSUME(size > 1);
+#if USE_AVL == 1
     /* FIXME: The AVL-tree is the bottleneck of this algorithm. A Treap
        [R. Seidel and C. R. Aragon. Randomized search trees.  Algorithmica,
        16:464–497, 1996] may be far more efficient by allowing to remove a
@@ -402,6 +430,42 @@ find_nondominated_3d_impl_sorted(const double ** restrict rows, size_t size,
                    : printf_point("insert before point: ", point, 3, "\n"));
             (++node)->item = pj;
             avl_insert_before(&tree, nodeaux, node);
+#else
+    TreapNode *tnodes = malloc(size * sizeof(*tnodes));
+    assert(tnodes != NULL);
+    TreapNode * node = tnodes;
+    double pk0 = rows[0][0], pk1 = rows[0][1], pk2 = rows[0][2];
+    treap_node_init(node, pk0, rows[0]);
+    Treap tree;
+    treap_init_with_single_node(&tree, node);
+    node++;
+
+    // In this context, size means "no dominated solution found".
+    size_t new_size = size;
+    bool prev_dominated = false;
+    for (size_t j = 1; j < size; ++j) {
+        const double * restrict pj = rows[j];
+        DEBUG2(printf_point("pj = [ ", pj, 3, " ], "));
+        const double pj0 = pj[0], pj1 = pj[1], pj2 = pj[2];
+        if ((pk0 > pj0) | (pk1 > pj1)) {
+            // Check if pj is dominated by a point in the tree.
+            /* In a valid 2-D frontier, x increases and y decreases.  Therefore
+               the only existing point that can dominate pj is the frontier
+               predecessor with the largest key <= pj0.  */
+            TreapNode *pred = treap_find_le(&tree, pj0);
+            if (pred != NULL && pred->item[1] <= pj1)
+                goto j_is_dominated;
+            /* pj is not dominated by an existing frontier point.
+
+               Insert it and detach every existing point dominated by it.
+
+               The returned treap contains exactly those displaced nodes, but
+               we do not need to traverse it here because those nodes will
+               never again participate in dominance queries.  */
+            treap_node_init(node, pj0, pj);
+            (void) treap_insert_and_displace(&tree, node);
+            node++;
+#endif
             // Fall-through to j_is_NOT_dominated.
         } // Handle duplicates and points that are dominated by the immediate previous one.
         else if (!keep_weakly // Don't keep duplicates.
@@ -430,7 +494,6 @@ find_nondominated_3d_impl_sorted(const double ** restrict rows, size_t size,
         rows[j] = NULL;
         new_size--;
     }
-
 early_end:
     free(tnodes);
     return new_size;