This Ocaml code computes the boundary of an abstract oriented simplicial complex. Such a complex partitions some space (manifold) into simplexes. Some finite set of points, which includes all the vertices of all the simplexes, is enumerated and this numbering is used to specify the vertices of each constituent simplex. Generally a list of integers (int list) here depicts some simplex by enumerating its vertices. The logic of this program has no notion of location for these points and the manifold is thus ‘abstract’.
We include in red the type assigned by ocaml to each of the program fragments that we have copied into this explanation. We use a sort routine that is expanded and expounded here. Here is a file with the boiled down code, shorn of prose.
exception Repeated_Vertex;; exception Logic_error;; let pe vl = match let rec mpa l th = match l with [] -> [], [] | (a, b)::c -> if th >= b then raise Repeated_Vertex; match mpa c b with sl, l -> (b::sl, (ref false, a)::l) in mpa (sort (fun x y -> match (x, y) with (_, b), (_, d) -> b < d) (let rec numr l n = match l with [] -> [] | a::b -> (n, a)::(numr b (n + 1)) in numr vl 0)) (-1) with svl, spl -> let ar = Array.of_list spl in let rec dp i p = if i < 0 then p else let next = dp (i-1) in match ar.(i) with a, n -> if !a then next (not p) else let rec lop v = match ar.(v) with w, x -> if !w then raise Logic_error; w := true; if v = i then () else lop x in lop n; next p in dp (Array.length ar - 1) false, svl;; pe : int list -> bool * int list pe [3; 6; 5; 8];; => - : bool * int list = (true, [3; 5; 6; 8])pe takes a simplex and sorts the vertex numbers and reports parity of the permutation of the sort. It is described in expanded form here.
let rec sl l = match l with [] -> [] | h::t -> t::List.map (fun q -> h::q) (sl t);; sl : 'a list -> 'a list list sl [3; 5; 6];; => - : int list list = [[5; 6]; [3; 6]; [3; 5]](sl z) produces a list of all sublists of z that lack just one member of z. If x is an n-simplex then (sl x) is a list of its n+1 faces, each of which is an (n−1)-simplex. The orientation of these faces alternates as we pass down this list. This accounts for some obscure code later.
let rc b = List.fold_left (fun acc fac -> match fac with (o, vl) -> (* Do all faces of this zone with face list fac. *) (* Yield a pair of face lists assorted by parity. *) match (List.fold_left (fun acc1 face -> match acc1 with od, (even, odd) -> (not od, if o then (odd, face::even) else (face::odd, even))) (false, acc) (sl vl)) with (o, (l, r)) -> if o then (r, l) else (l, r)) ([], []) (List.map pe b);; rc : int list list -> int list list * int list list rc [[1; 2; 3; 4]; [5; 3; 2; 4]];; => - : int list list * int list list = ([[2; 3; 4]; [2; 4; 5]; [1; 2; 3]; [1; 3; 4]], [[2; 3; 5]; [3; 4; 5]; [1; 2; 4]; [2; 3; 4]])rc takes a list of simplexes and considers each face of each simplex. It accumulates these in two lists assorted by parity. Notice logic to alternate parity of faces of one simplex. Note in the given example of two tetrahedra that share the face [2; 3; 4], that the computed result includes that shared face in each list.
let b1 q = match (rc q) with a, b -> let ob = List.sort compare in (ob a, ob b);; b1 : int list list -> int list list * int list list b1 [[1; 2; 3; 4]; [5; 3; 2; 4]];; => - : int list list * int list list = ([[1; 2; 3]; [1; 3; 4]; [2; 3; 4]; [2; 4; 5]], [[1; 2; 4]; [2; 3; 4]; [2; 3; 5]; [3; 4; 5]])b1 applies rc to a set of simplexes and sorts the results, better to find matches between them. List.sort compare produces a sort for simplexes.
let rec elim q = match q with | [], [] -> [], [] | a, [] -> a, [] | [], b -> [], b | (a::b as rev), (c::d as obv) -> let w = compare a c in if w < 0 then match elim (b, obv) with p, q -> a::p, q else if w > 0 then match elim (rev, d) with p, q -> p, c::q else elim (b, d);; elim : 'a list * 'a list -> 'a list * 'a list elim ([[1; 2; 3]; [1; 3; 4]; [2; 3; 4]; [2; 4; 5]], [[1; 2; 4]; [2; 3; 4]; [2; 3; 5]; [3; 4; 5]]);; =>- : int list list * int list list = ([[1; 2; 3]; [1; 3; 4]; [2; 4; 5]], [[1; 2; 4]; [2; 3; 5]; [3; 4; 5]])elim takes a pair of sorted simplex lists, and eliminates those that are in both. Notice that the shared face [2; 3; 4] has been eliminated in the result. Such simplexes are faces between neighboring zones which are properly oriented so as not to belong to the boundary of the complex.
let comb a = elim (b1 a);; comb : int list list -> int list list * int list listcomb merely combines the sorting and elimination operations.
let boundary x = let rec turn (p, q) = match q with [] -> p | r::s -> turn (((match r with a::b::c -> b::a::c)::p), s) in turn (comb x);; boundary : int list list -> int list list boundary [[1; 2; 3; 4]; [5; 3; 2; 4]];; => - : int list list = [[4; 3; 5]; [3; 2; 5]; [2; 1; 4]; [1; 2; 3]; [1; 3; 4]; [2; 4; 5]]boundary combines the previous processing and turns the pair of lists into one list, by flipping a pair of vertices in each member of the first list, thus making the output format into the same as the input. We can now illustrate the famous theorem ‘the boundary of the boundary is 0’:
boundary (boundary [[1; 2; 3; 4]; [5; 3; 2; 4]]);; =>- : int list list = []If you want to learn of repeated boundary members in a sorted list of simplexes try:
let rec rep x = match x with [] -> [] | [a] -> [] | a::b::c -> let z = rep (b::c) in if a = b then a::z else z;; rep : 'a list -> 'a list