This code computes moments 0 and 1 of general polyhedra but the topology of the surface must be done manually.
The routine <tt>boundary</tt> is merely a good check on the topology.
We have thus another way of computing moments corresponding to c((1, u, v)) where 0 ≤ u, 0 ≤ v and u+v ≤ 1.
The clipped cube still has 8 vertices and unchanged surface topology.
Numbering the vertices conventionally the OCaml list <tt>cpm</tt> below enumerates the triangles that together bound the clipped cube.
This provides the topological relations between the simplexes that bound the body.
<p>With all of <a href=../IP/ocaml/bndry.ml>boundary logic</a> and <a href=../../code/ocaml/Linear.ml>Matrix stuff</a>
<pre>let cmp = [[0; 1; 2]; [3; 2; 1]; [4; 6; 5]; [7; 5; 6];
[1; 0; 4]; [1; 4; 5]; [0; 2; 6]; [0; 6; 4];
[2; 3; 6]; [6; 3; 7]; [7; 3; 1]; [7; 1; 5]];;</pre>
Observe that <tt>boundary cmp;; ⇒ int list list []</tt>.
This is good evidence that we got the triangle list right.
Array <tt>geo</tt> is the coordinates of the vertices of the clipped cube.

<pre>let geo = let u = 0.39 and v = 0.47 in
 [|[|0.; 0.; 0.|]; [|0.; 0.; 1.|]; [|0.; 1.; 0.|]; [|0.; 1.; 1.|];
  [|1.; 0.; 0.|]; [|1. -. v; 0.; 1.|];
  [|1. -. u; 1.; 0.|]; [|1. -. u -. v; 1.; 1.|]|];;

let ai = Array.init in
let mq f = let wtt = Array.of_list f in
  ai (List.length f) (fun j -> geo.(wtt.(j))) in
let rec fac n = if n=0 then 1. else (float n) *. (fac (n-1)) in
let n = (List.length (List.hd cmp)) in
let fn = fac n and (m0, m1) = let az = Array.make n 0. in
List.fold_left (fun (cm0, cm1) zn ->
  let matr = mq zn in
  let ax = (Array.fold_left Linear.vadd az matr)
  and det = Linear.determ matr in
        cm0 +. det, ai n (fun j -> det *. ax.(j) +. cm1.(j)))
    (0., az) cmp in
(m0 /. fn, Linear.sm (1. /. ((float (n+1)) *. fn)) m1);;
</pre>
<hr><pre>cube corners: (0, 0, 0), (0, 0, 1), (0, 1, 0), (0, 1, 1),
(1, 0, 0), (1-v, 0, 1), (1-u, 1, 0), (1-u-v, 1, 1)</pre>