Initial commit

.gitignore

+_build/

bin/bench/dune

+(executable
+ (public_name robench)
+ (name robench)
+ (libraries 
+   DataStruct 
+   CopointLib 
+   DissimilarityLib 
+   RandomDissimilarity 
+   Recognition
+   core core_unix.command_unix core_bench
+ )
+)

bin/bench/robench.ml

+open! Core_bench
+open! CopointLib
+open! RandomDissimilarity
+open! DissimilarityLib
+
+
+let imperative_copoint_algorithm diss =
+  diss
+  |> ImperativeCCNP.find_compatible_order
+  |> ignore
+
+
+let functional_copoint_algorithm diss =
+  diss 
+  |> RobinsonCCNP.find_compatible_order
+  |> ignore
+
+let pivotpair_algorithm diss = 
+  Recognition.Pivotpair.algo diss (fun i -> i)
+  |> ignore
+
+
+let toeplitz ~k dim = 
+  Toeplitz.toeplitz012 ~n:dim ~k:(min k (dim-2))
+    
+let incrementing dim = 
+  Incrementing.random (Random.get_state ()) dim
+
+let from_pqtree dim = 
+  FromPqtree.uniform dim 
+
+let from_intervals ~k dim = 
+  FromPqtree.from_intervals ~n:dim ~k
+
+let from_rectangle ~height ~width dim = 
+  FromModel.from_points_in_rectangle ~width ~height ~count:dim
+
+let dimensions = [10;100;1000]
+
+type 'a check = Y of 'a  | N of 'a 
+
+let _ = N 0 
+
+let algorithms = 
+  [
+    Y ("copoint-imp", imperative_copoint_algorithm);
+    Y ("copoint-fun", functional_copoint_algorithm);
+    Y ("pivotpair", pivotpair_algorithm)
+  ]
+
+
+let generators =
+  [
+    Y ("toeplitz1", toeplitz ~k:1);
+    Y ("toeplitz2", toeplitz ~k:2);
+    Y ("toeplitz5", toeplitz ~k:5);
+    Y ("incrementing", incrementing);
+    Y ("fromPQtree", from_pqtree);
+    Y ("fromIntervals5", from_intervals ~k:5);
+    Y ("fromIntervals10", from_intervals ~k:10);
+    Y ("fromRectangle100", from_rectangle ~width:100. ~height:10.);
+    Y ("fromRectangle1000", from_rectangle ~width:1000. ~height:10.);
+    Y ("fromRectangle10000", from_rectangle ~width:10000. ~height:10.);
+  ]
+
+let select check_list =
+  List.filter_map (function Y a -> Some a | N _ -> None) check_list 
+
+let tests = 
+  let (>>=) e f = DataStruct.MoreList.flat_map f e in 
+  select generators >>= fun (gen_name, gen) -> 
+  select algorithms >>= fun (algo_name, algo) -> 
+  [ 
+    Bench.Test.create_indexed 
+      ~name:(Format.sprintf "%s %s" algo_name gen_name) 
+      ~args:dimensions 
+      (fun dim ->
+        Core.Staged.stage 
+          (fun () -> dim |> gen |> Dissimilarity.shuffle |> algo)
+      )
+  ]
+
+let () = 
+  tests
+  |> Bench.make_command
+  |> Command_unix.run
+       

bin/dune

+(executable
+ (public_name copoint)
+ (name main)
+ (libraries RobinsonLib CopointLib DissimilarityLib Recognition))

bin/main.ml

+open DissimilarityLib
+open CopointLib
+
+(* let _list = 
+ *   [ 0; 3; 4; 1;
+ *     3; 0; 1; 2;
+ *     4; 1; 0; 3;
+ *     1; 2; 3; 0 
+ *   ]
+ * 
+ * 
+ * let print_ordered_matrix (module D : Dissimilarity.D with type t = int) order = 
+ *   List.iter (fun i -> 
+ *       List.iter (fun j -> 
+ *           Printf.printf "%d " (D.d i j)
+ *         )
+ *         order;
+ *       Printf.printf "\n"
+ *     )
+ *     order *)
+
+
+
+
+let _from_stdin algo =  
+  match IO.Plain.from_stdin () with 
+  | None -> Printf.printf "Parsing error"
+  | Some diss ->
+     let order = algo diss in
+     let pp_sep fmt () =  Format.pp_print_string fmt ", " in
+     Format.printf "Order : %a\n%!" 
+       Format.(pp_print_option (pp_print_list ~pp_sep pp_print_int))
+       order;
+     Option.iter (fun elements -> 
+         Format.printf "Matrix:@,@[<v 2>%a@]\n%!"
+           IO.Plain.pretty_printer 
+           Dissimilarity.({ diss with elements })
+       )
+       order;  
+     Format.printf "%!";
+     exit 0
+
+
+
+let _main = 
+  let _algo1 diss = RobinsonCCNP.find_compatible_order diss in 
+  let _algo2 diss = ImperativeCCNP.find_compatible_order diss |> Option.map Array.to_list in 
+  let _algo3 diss = Recognition.Pivotpair.algo diss (fun i -> i) in
+  _from_stdin _algo1 

dune-project

+(lang dune 3.0)
+
+(name robinson)
+
+(generate_opam_files true)
+
+(source
+ (github username/reponame))
+
+(authors "Guyslain Naves")
+
+(maintainers "Guyslain Naves")
+
+(license LICENSE)
+
+(documentation https://url/to/documentation)
+
+(package
+ (name robinson)
+ (synopsis "Algorithms for Robinson Dissimilarities")
+ (description "A longer description")
+ (depends ocaml dune fix feat RobinsonLib)
+ (tags
+  (topics "to describe" your project)))
+
+; See the complete stanza docs at https://dune.readthedocs.io/en/stable/dune-files.html#dune-project

lib/combi/combinatorics.ml

+open! DataStruct
+
+let one = Big_int_Z.big_int_of_int 1 
+
+module RND = FeatCore.RandomBigInt.Make(Feat.Num)(Random : FeatCore.RandomSig.S)
+
+let factorials = 
+  let open LazyList.Infinite in 
+  let open Big_int_Z in 
+  unfold (fun (n,f) -> (n+1, mult_big_int (big_int_of_int (n+1)) f)) (1,one) 
+  |> map snd 
+
+let factorial n = 
+  LazyList.Infinite.get n factorials
+
+
+let rec random_permutation = function
+  | [] -> []
+  | [ _ ] as l -> l 
+  | list -> 
+     let (left,right) = List.partition (fun _ -> Random.bool ()) list in
+     List.rev_append (random_permutation left) (random_permutation right)
+
+
+let choose choices = 
+  let sum = 
+    choices
+    |> List.map snd 
+    |>  List.fold_left Feat.Num.add Feat.Num.zero
+  in 
+  let r = RND.random sum in 
+  let rec go r = function
+    | (hd,k)::_ when r < k -> hd 
+    | (_,k)::tail ->  go (Feat.Num.sub r k) tail 
+    | [] -> assert false
+  in  
+  go r choices
+
+
+let random_sublist list = 
+  let l = List.length list in 
+  let a = Random.int l in 
+  let b = Random.int l in 
+  let (mini,maxi) = (min a b, max a b) in 
+  list
+  |> MoreList.drop mini
+  |> MoreList.take (maxi - mini + 1)
+
+(* module GenSeries = 
+ *   struct 
+ * 
+ *     type t =  
+ *       int LazyList.Infinite.t
+ *     
+ *     open LazyList.Infinite
+ *     open Big_int_Z 
+ *        
+ *     let (++) = Big_int_Z.add_big_int
+ *     let ( ** ) = Big_int_Z.mult_big_int   
+ * 
+ * 
+ *     let zero = unfold (fun _ -> zero_big_int) zero_big_int
+ *     let one = unfold (fun _ -> zero_big_int) (big_int_of_int 1)
+ *     let unit = lazy (Cons (zero_big_int, one))
+ * 
+ *     let from_general_form f =
+ *       map f ints 
+ * 
+ *     let increment a = lazy (Cons (zero_big_int, a))
+ * 
+ *     let sum a b = zip_with (++) a b 
+ * 
+ *     let product a b = 
+ *       ints 
+ *       |> map (fun n -> 
+ *              List.combine
+ *                (a |> take (n + 1))
+ *                (List.rev (b |> take (n + 1)))
+ *              |> List.map (fun (a,b) -> a ** b) 
+ *              |> List.fold_left (++) zero_big_int 
+ *            )
+ * 
+ * 
+ *     let list a = 
+ *       unfold (fun (previous, _) ->
+ *           let next = 
+ *             a 
+ *             |> drop 1
+ *             |> take (List.length previous)
+ *             |> List.combine previous
+ *             |> List.map (fun (pi, aj) -> pi ** aj)
+ *             |> List.fold_left (++) zero_big_int
+ *           in 
+ *           (next::previous, next)
+ *         ) 
+ *         ( [big_int_of_int 1], big_int_of_int 1)
+ *       |> map snd 
+ * 
+ * 
+ *     let non_empty_list a =
+ *       product a (list a)
+ * 
+ *     let rec list_lb lb a =
+ *       if lb = 0 then list a 
+ *       else product a (list_lb (lb-1) a)
+ * 
+ * 
+ *   end *)

lib/combi/combinatorics.mli

+val factorial : int -> Big_int_Z.big_int
+
+
+val random_permutation : 'elt list -> 'elt list
+
+val random_sublist : 'elt list -> 'elt list 
+
+
+val choose : ('elt * Big_int_Z.big_int) list -> 'elt 

lib/combi/dune

+(library
+ (name Combi)
+ (libraries feat fix DataStruct)
+)

lib/combi/enumeration.ml

+
+let int_enumerator start =
+  Seq.unfold (fun b -> Some (b, b+1)) start
+  
+
+let rec valuation_2 n = 
+  if n mod 2 = 0 then 1 + valuation_2 (n/2) 
+  else 0
+
+
+let gray_code_enumerator = 
+  Seq.map valuation_2 (int_enumerator 1)
+
+
+let subset_enumerator elements =
+  let n = List.length elements in 
+  let rec flip k elements subset = 
+    match elements, subset with  
+    | e1::_, s::sub when k = 0 && e1 = s -> sub 
+    | e1::_, _ when k = 0 -> e1::subset 
+    | e1::elts, s::sub when e1 = s -> s :: flip (k-1) elts sub 
+    | _::elts, _ -> flip (k-1) elts subset
+    | _ -> assert false 
+  in 
+  let step (index,subset) =
+        let two_val = valuation_2 index in 
+        if two_val >= n then None
+        else 
+          let next_subset = flip two_val elements subset in 
+          Some (next_subset, (index+1, next_subset))
+  in
+  Seq.cons elements (Seq.unfold step (1, []))
+  
+
+
+let rec fact_valuation f n =
+  let r = n mod f in 
+  if r = 0 then fact_valuation (f+1) (n/f) 
+  else (f,r)
+
+let transposition k = 
+  let (n,r) = fact_valuation 2 k in 
+  if n mod 2 = 0 then (r,n) 
+  else (1,n) 
+
+let heap_transposition_enumerator = 
+  Seq.map transposition (int_enumerator 1)
+
+
+let permutation_enumerator elements = 
+  let n = List.length elements in 
+  let rec replace i e = function
+    | f::tail when i = 0 -> (f,e::tail)
+    | head::tail -> 
+       let (f,tail') =  replace (i-1) e tail in 
+       (f, head::tail')
+    | [] -> assert false
+  in 
+  let rec transpose i j = function 
+    | ei::tail when i = 0 -> 
+       let (ej, tail') = replace (j-1) ei tail in 
+       ej::tail'
+    | head::tail -> 
+       head :: transpose (i-1) (j-1) tail 
+    | [] -> assert false
+  in 
+  let step (index, permutation) =
+    let (i,j) = transposition index in 
+    if j > n then None
+    else 
+      let next_permutation = transpose (i-1) (j-1) permutation in 
+      Some (next_permutation, (index+1, next_permutation))
+  in 
+  Seq.cons elements (Seq.unfold step (1, elements))
+  
+
+
+let pair_enumerator enum_a enum_b = 
+  let open Seq in 
+  let rec step (state_a, current_b, state_b) = 
+    match state_a () with 
+    | Nil -> 
+       begin match state_b () with 
+       | Nil -> None 
+       | Cons (new_b, new_state_b) -> step (enum_a, new_b, new_state_b)
+       end 
+    | Cons (new_a, new_state_a) -> 
+       Some ((new_a,current_b), (new_state_a, current_b, state_b))
+  in 
+  match enum_a (), enum_b () with 
+  | Cons (a, state_a), Cons (b, state_b) -> 
+     Seq.cons (a,b) (Seq.unfold step (state_a,b,state_b))
+  | _,_ -> Seq.empty
+
+
+let (>>=) enum_a a_to_enum_b =
+  Seq.flat_map a_to_enum_b enum_a
+  
+let rec list_enumerator = function 
+  | [] -> Seq.return []
+  | enum_head :: others ->
+     list_enumerator others >>= fun tail ->
+     enum_head >>= fun head -> 
+     Seq.return (head::tail)
+
+
+let permutation_list_enumerator enum_list =
+  permutation_enumerator enum_list >>= list_enumerator 
+  
+
+
+(* let rec take n seq = 
+ *   match seq () with 
+ *   | _ when n = 0 -> []          
+ *   | Seq.Nil -> []
+ *   | Seq.Cons (e,seq) -> e :: take (n-1) seq 
+ * 
+ * 
+ * let _ = 
+ *   subset_enumerator [1;2;3;4;5;6] |> List.of_seq |> List.length
+ * 
+ * let _ = 
+ *   permutation_enumerator [1;2;3;4;5] |> List.of_seq |> List.length
+ * 
+ * let _ =
+ *   pair_enumerator (List.to_seq [1;2;3;4]) (List.to_seq ['a';'b';'c']) |> List.of_seq
+ * 
+ * let _ = 
+ *   list_enumerator [ List.to_seq [1;2;3]; List.to_seq [4;5]; List.to_seq [6]; List.to_seq [7;8] ] |> List.of_seq
+ * 
+ * let _ = 
+ *   permutation_list_enumerator [ List.to_seq [1;2;3]; List.to_seq [4;5]; List.to_seq [6]; List.to_seq [7;8] ] |> List.of_seq *)

lib/combi/enumeration.mli

+val int_enumerator : int -> int Seq.t 
+
+val gray_code_enumerator : int Seq.t 
+
+val subset_enumerator : 'elt list -> 'elt list Seq.t 
+
+val heap_transposition_enumerator : (int * int) Seq.t
+
+val permutation_enumerator : 'elt list -> 'elt list Seq.t 
+
+
+val pair_enumerator : 'a Seq.t -> 'b Seq.t -> ('a * 'b) Seq.t 
+
+val list_enumerator : 'e Seq.t list -> 'e list Seq.t 
+
+val (>>=) : 'a Seq.t -> ('a -> 'b Seq.t) -> 'b Seq.t 
+
+val permutation_list_enumerator : 'e Seq.t list -> 'e list Seq.t 
+

lib/copoint/dune

+(library
+ (name CopointLib)
+ (libraries DissimilarityLib)
+)

lib/copoint/imperativeCCNP.ml

+open! DataStruct 
+open DissimilarityLib
+
+module IMap = Map.Make(Int) 
+
+
+
+module Make = 
+   functor (D : sig type t val d : t -> t -> int end) ->
+   struct 
+
+     let refine q set =
+       let tree = ref IMap.empty in 
+       for i = 0 to Array.length set - 1 do 
+         let dist = D.d q set.(i) in 
+         if IMap.mem dist !tree then 
+           tree := IMap.add dist (set.(i) :: IMap.find dist !tree) !tree 
+         else 
+           tree := IMap.add dist [set.(i)] !tree 
+       done;
+       IMap.to_seq !tree 
+       |> Seq.map snd 
+       |> Array.of_seq
+
+
+     let swap arr i j =
+       let old_i = arr.(i) in 
+       arr.(i) <- arr.(j);
+       arr.(j) <- old_i
+
+
+     let reverse_short_distance p q subsets = 
+       let is_close subset = 
+         subset <> [] && D.d q (List.hd subset) <= D.d p q
+       in 
+       let i = ref 0 in 
+       while !i < Array.length subsets && is_close subsets.(!i) do incr i done;
+       for j = 0 to !i/2 - 1 do 
+         swap subsets j (!i-j-1)
+       done
+
+
+     let choose_pivot inners outers = 
+       if inners = [] then 
+         (List.hd outers, inners, List.tl outers, true)
+       else 
+         (List.hd inners, List.tl inners, outers, false)
+
+
+     let rec recursive_refine p (inners, set, outers) =
+       if (inners = [] && outers = []) || Array.length set <= 1 then 
+         [set]
+       else 
+         let (q, ins, outs, q_is_in_outers) = choose_pivot inners outers in 
+         let subsets = refine q set in 
+         if q_is_in_outers then reverse_short_distance p q subsets;
+         let concats = ref [] in 
+         for i = Array.length subsets - 1 downto 0 do 
+           let ins_i, outs_i = ref ins, ref outs in 
+           for j = 0 to i-1 do 
+             ins_i := List.rev_append subsets.(j) !ins_i;
+           done;
+           for j = i+1 to Array.length subsets - 1 do 
+             outs_i := List.rev_append subsets.(j) !outs_i;
+           done;
+           let copoints = recursive_refine p (!ins_i, Array.of_list subsets.(i),!outs_i) in 
+           concats := List.append copoints !concats;
+         done;
+         !concats
+
+     exception Not_Robinson 
+
+     
+     let separate_if_separable p copoint =
+       let l = Array.length copoint in  
+       let x_min = copoint.(0) in 
+       let x_max = copoint.(l-1) in 
+       let delta = D.d p x_min in 
+       let diam = D.d x_min x_max in 
+       if diam <= delta then 
+         [(x_min, copoint)]
+       else 
+         let i = ref 0 in 
+         while not (
+                   D.d x_min copoint.(!i) <= delta  
+                   && D.d x_max copoint.(!i+1) <= delta 
+                   && D.d copoint.(!i) copoint.(!i+1) >= delta
+                 )
+         do 
+           incr i;
+           if !i >= l - 1 then raise Not_Robinson
+         done;
+         [ (x_min, Array.sub copoint 0 (!i + 1));
+           (x_max, Array.sub copoint (!i + 1) (l - !i - 1)) 
+         ]
+               
+
+     type side = Left | Right 
+
+     let append_all array list_ref =
+       for i = 0 to Array.length array - 1 do 
+         list_ref := array.(i) :: !list_ref
+       done
+       
+
+     let sort_by_bipartition p set =
+       let l = Array.length set in 
+       let sides = Array.make l Right in 
+       let lefts = ref [] in 
+       let rights = ref [] in 
+       let last_undecided = ref (l-1) in 
+       let move_right i = 
+         append_all (snd set.(i)) rights;
+         sides.(i) <- Right;
+         decr last_undecided
+       in 
+       let move_left i = 
+         append_all (snd set.(i)) lefts;
+         sides.(i) <- Left;
+         decr last_undecided
+       in
+       for i = l-1 downto 0 do 
+         let q = fst set.(i) in 
+         if i = !last_undecided then move_right i;
+         for j = !last_undecided downto 0 do 
+           let x = fst set.(j) in 
+           if not (D.d q x = D.d p q) then 
+             if (D.d q x < D.d p q) = (sides.(i) = Right) then begin
+                 for k = !last_undecided downto j+1 do move_left k done;
+                 move_right j;
+               end
+             else begin
+                 for k = !last_undecided downto j+1 do move_right k done;
+                 move_left j
+               end;
+         done;
+       done;
+       List.rev_append !lefts (p :: !rights)
+
+
+     let rec find_compatible_order set =
+       let n = Array.length set in 
+       if n <= 1 then set 
+       else 
+         let p = set.(0) in 
+         let set' = Array.sub set 1 (n-1) in  
+         let copoints = recursive_refine p ([p], set', []) |> Array.of_list in 
+         let represented_copoints = ref [] in 
+         for i = 0 to Array.length copoints - 1 do 
+           let order_i = find_compatible_order copoints.(i) in 
+           represented_copoints := List.rev_append (separate_if_separable p order_i) !represented_copoints;
+         done;
+         sort_by_bipartition p (Array.of_list (List.rev !represented_copoints))
+         |> Array.of_list
+   end
+
+
+let find_compatible_order (type elt) diss = 
+  let open Dissimilarity in 
+  let module D = struct type t = elt let d = diss.d end in 
+  let module M = Make(D) in 
+  try 
+    Some (M.find_compatible_order (Array.of_list diss.elements))
+  with 
+    | M.Not_Robinson -> None 

lib/copoint/imperativeCCNP.mli

+open DissimilarityLib 
+
+
+val find_compatible_order : 'elt Dissimilarity.t -> 'elt array option
+
+
+

lib/copoint/robinsonCCNP.ml

+open! DataStruct
+open DissimilarityLib 
+
+module IMap = Map.Make(Int)
+
+module Queue = FQueue 
+
+
+
+type side = Right | Left
+
+
+
+module Make = 
+  functor (D: sig type t val d : t -> t -> int end) -> 
+  struct 
+    
+
+    let refine q set = 
+      let insert tree element =
+        let dist = D.d q element in 
+        let previous_list = IMap.find_opt dist tree |> Option.value ~default:[] in 
+        IMap.add dist (List.cons element previous_list) tree
+      in  
+      set 
+      |> List.fold_left insert IMap.empty
+      |> IMap.to_seq
+      |> Seq.map snd
+      |> List.of_seq
+
+
+    let reverse_short_distance p q subsets = 
+      let is_close = function
+        | x::_ -> D.d q x <= D.d p q
+        | [] -> false
+      in
+      let (close_sets,far_sets) =
+        MoreList.take_while ~f:is_close subsets 
+      in
+      List.rev_append close_sets far_sets 
+        
+
+      
+    let rec recursive_refine p (inners, set, outers) =
+      match inners, outers with
+      | [], [] -> [set]
+      | _ when List.length set <= 1 -> [set]
+      | (q::ins,outs) 
+      | ([] as ins, q::outs) -> 
+         let subsets = refine q set in 
+         let subsets' = 
+           if inners = [] then 
+             reverse_short_distance p q subsets 
+           else 
+             subsets 
+         in
+         let recurse (subset, previous_subsets, next_subsets) =
+           recursive_refine p 
+             (List.rev_append (List.flatten previous_subsets) ins,
+              subset,
+              List.rev_append (List.flatten next_subsets) outs
+             ) 
+         in
+         subsets' 
+         |> MoreList.extend
+         |> List.map recurse  
+         |> List.flatten 
+         
+         
+    exception Not_Robinson
+
+
+    let separate_if_separable p = function
+      | [] -> []
+      | (x_min::_) as copoint ->  
+         let x_max = List.fold_left (fun _ last -> last) x_min copoint in 
+         let diameter = D.d x_min x_max in 
+         let delta = D.d x_min p in 
+         if diameter <= delta then [(x_min,copoint)]
+         else 
+           let rec split prefix = function
+             | x::y::suffix when D.d x_min x <= delta && D.d y x_max <= delta && D.d x y >= delta ->
+                [ (x_min, x::prefix); (x_max, y::suffix) ]
+             | x::suffix -> split (x::prefix) suffix
+             | [] -> raise Not_Robinson
+           in 
+           split [] copoint
+
+
+
+    type representant = D.t * D.t list
+    type dispatch = {
+        lefts : representant list;
+        rights : representant list;
+        pivots : (representant * side) Queue.t; 
+        skipped : representant list
+      }
+    let insert_to_left element dispatch =
+      { dispatch with 
+       lefts = element::dispatch.lefts;
+       pivots = Queue.enqueue dispatch.pivots (element,Left)
+      }
+    let insert_to_right element dispatch =
+      { dispatch with
+        rights = element::dispatch.rights;
+        pivots = Queue.enqueue dispatch.pivots (element,Right)
+      }
+    let skip element dispatch = 
+      { dispatch with skipped = element :: dispatch.skipped }
+
+    let move_skipped_to_left dispatch = 
+      { dispatch with  skipped = [] }
+      |> List.fold_right insert_to_left dispatch.skipped
+    let move_skipped_to_right dispatch =
+      { dispatch with skipped = [] }
+      |> List.fold_right insert_to_right dispatch.skipped 
+    let reverse_skipped dispatch =
+      { dispatch with skipped = List.rev dispatch.skipped }
+
+    let get_next_pivot dispatch = 
+      match Queue.view dispatch.pivots, dispatch.skipped with
+      | Some ((q,side), pivots),_ -> (q, side, { dispatch with pivots })
+      | None, q::skipped -> (q, Right, insert_to_right q { dispatch with skipped })
+      | None, [] -> assert false
+      
+    let rec repeat_until ~condition ~f start = 
+      if condition start then start
+      else repeat_until ~condition ~f (f start)
+
+    let sort_by_bipartition p set = 
+      let bipart_using_pivot (q, q_side, dispatch) = 
+        List.fold_left (fun dispatch x ->
+            if D.d (fst x) (fst q) = D.d p (fst q) then skip x dispatch
+            else if (D.d (fst x) (fst q) < D.d p (fst q)) = (q_side = Right) then 
+              dispatch |> move_skipped_to_left |> insert_to_right x
+            else 
+              dispatch |> move_skipped_to_right |> insert_to_left x
+          )
+          { dispatch with skipped = [] }
+          dispatch.skipped
+      in
+      let { lefts; rights; _ } =
+        repeat_until
+          ~condition:(fun dispatch -> dispatch.skipped = [])
+          ~f:(fun dispatch ->
+            get_next_pivot dispatch
+            |> bipart_using_pivot 
+            |> reverse_skipped            
+          )
+          { lefts = []; rights = []; pivots = Queue.empty; skipped = List.rev set }
+      in
+      List.rev_append (List.map snd lefts) ([p] :: List.map snd rights)
+      |> List.flatten
+
+
+      let rec find_compatible_order = function
+        | [] -> []
+        | p::set -> 
+           let copoints = recursive_refine p ([p],set,[]) in 
+           copoints 
+           |> List.map find_compatible_order
+           |> List.concat_map (separate_if_separable p)
+           |> sort_by_bipartition p 
+           
+  end
+
+
+let find_compatible_order (type elt) diss = 
+  let open Dissimilarity in 
+  let module D = 
+    struct 
+      type t = elt 
+      let d = diss.d 
+    end 
+  in 
+  let module M = Make(D) in 
+  try
+    Some (M.find_compatible_order diss.elements)
+  with 
+  | M.Not_Robinson -> None

lib/copoint/robinsonCCNP.mli

+open DissimilarityLib 
+
+val find_compatible_order : 'elt Dissimilarity.t -> 'elt list option

lib/datastruct/dune

+(library
+ (name DataStruct)
+)

lib/datastruct/fQueue.ml

+
+type 'elt t = 'elt list * 'elt list 
+
+let empty = ([], [])
+
+let enqueue queue elt = match queue with
+  | ([],[]) -> ([elt],[]) 
+  | (pre,post) -> (pre,elt::post) 
+  
+let is_empty (pre,_) = pre = [] 
+
+let view = function 
+  | ([],_) -> None 
+  | ([first],others) -> Some (first, (List.rev others,[]))
+  | (first::pre,post) -> Some (first, (pre,post))
+
+let of_list list = (List.rev list, [])

lib/datastruct/fQueue.mli

+type 'elt t 
+
+val empty : 'elt t 
+
+val enqueue : 'elt t -> 'elt -> 'elt t 
+
+val view : 'elt t -> ('elt * 'elt t) option 
+
+val is_empty : 'elt t -> bool 
+
+val of_list : 'elt list -> 'elt t 

lib/datastruct/lazyList.ml

+module PossiblyFinite =
+  struct 
+    
+    type 'elt t = 'elt node Lazy.t 
+    and 'elt node = 
+      | Nil 
+      | Cons of 'elt * 'elt t 
+
+
+    let rec unfold f init = 
+      lazy (
+          match f init with 
+          | None -> Nil 
+          | Some elt -> Cons (elt, unfold f elt)
+        )
+      
+    let view (lazy v) = v 
+
+    let rec map f llist =
+      lazy (
+          match view llist with 
+          | Nil -> Nil 
+          | Cons (elt, tail) -> Cons (f elt, map f tail)
+        )                
+
+    let rec zip_with f llist1 llist2 =
+      lazy (
+          match view llist1, view llist2 with 
+          | Nil, _ | _, Nil -> Nil 
+          | Cons (elt1,tail1), Cons (elt2, tail2) -> Cons (f elt1 elt2, zip_with f tail1 tail2)
+        )
+        
+
+    let rec take n llist = 
+      if n = 0 then []
+      else match view llist with 
+           | Nil -> [] 
+           | Cons (head,tail) -> head :: take (n-1) tail 
+
+
+    let rec drop n llist =
+      if n <= 0 then llist 
+      else 
+        match view llist with 
+        | Nil -> llist 
+        | Cons (_,tail) -> drop (n-1) tail 
+
+    let head llist = 
+      match view llist with 
+      | Nil -> None 
+      | Cons (head, _) -> Some head 
+
+
+    let get n llist = 
+      llist 
+      |> drop n
+      |> head
+
+
+    let rec append list llist = 
+      match list with 
+      | [] -> llist 
+      | head :: tail -> lazy (Cons (head, append tail llist))
+
+
+    let rec rev_append list llist = 
+      match list with 
+      | [] -> llist 
+      | head::tail -> lazy (Lazy.force (rev_append tail (lazy (Cons (head, llist)))))
+
+
+    let to_seq llist = 
+      Seq.unfold (fun llist -> 
+          match view llist with 
+          | Nil -> None 
+          | Cons (head,tail) -> Some (head,tail)
+        ) llist 
+
+
+    let rec of_seq seq =
+      lazy (match seq () with 
+            | Seq.Nil -> Nil 
+            | Seq.Cons (head,tail) -> Cons (head, of_seq tail)
+        )
+
+    let range mini maxi = 
+      unfold (fun i -> if i >= maxi then None else Some (i+1)) (mini-1) 
+
+    let ints = 
+      unfold (fun i -> Some (i+1)) (-1)
+
+  end
+
+module Infinite = 
+  struct 
+
+    type 'elt t = 'elt node Lazy.t 
+    and 'elt node = Cons of 'elt * 'elt t 
+
+    let rec unfold f init = 
+      lazy (Cons (init, unfold f (f init)))
+
+    let rec unfoldl f linit = 
+      lazy (
+          let lazy init = linit in 
+          Cons (init, unfoldl f (lazy (f init)))
+        )
+
+    let rec const c = 
+      lazy (Cons (c, const c))
+
+    let rec append list llist =
+      match list with 
+      | [] -> llist 
+      | head::tail -> lazy (Cons (head, append tail llist))
+
+    let view (lazy list) = list 
+
+    let rec map f llist = 
+      lazy (
+          let Cons (head,tail) = view llist in 
+          Cons (f head, map f tail)
+        )
+
+    let rec zip_with f llist1 llist2 = 
+      lazy (
+          let lazy (Cons (elt1, tail1)) = llist1 in 
+          let lazy (Cons (elt2, tail2)) = llist2 in 
+          Cons (f elt1 elt2, zip_with f tail1 tail2)
+        )
+
+    let rec drop n llist = 
+      if n <= 0 then llist 
+      else 
+        lazy (
+            let Cons (_,tail) = view llist in 
+            view (drop (n-1) tail)
+          )
+
+    let rec take n llist = 
+      if n <= 0 then []
+      else 
+        let Cons (head,tail) = view llist in 
+        let tail' = take (n-1) tail in 
+        head::tail' 
+
+    let rec lazy_take n llist = 
+      lazy  (
+          if n <= 0 then PossiblyFinite.Nil
+          else 
+            let lazy (Cons (head,tail)) = llist in 
+            PossiblyFinite.Cons (head, lazy_take (n-1) tail)
+        )
+
+    let head (lazy (Cons (head,_))) = head 
+
+    let get n llist = head (drop n llist) 
+
+    let rec to_seq llist () = 
+      let lazy (Cons (head,tail)) = llist in 
+      Seq.Cons (head, to_seq tail)
+
+    let ints = unfoldl (fun n -> n + 1) (lazy 0) 
+
+    let rec fix f = f (fix f)
+
+  end