diff -r 495448cca603 -r 6aa955ad8122 src/sets/cube.rs
--- a/src/sets/cube.rs	Thu May 01 08:40:33 2025 -0500
+++ b/src/sets/cube.rs	Thu May 01 13:06:58 2025 -0500
@@ -16,25 +16,19 @@
 ```
 */
 
-use serde::ser::{Serialize, Serializer, SerializeTupleStruct};
-use crate::types::*;
 use crate::loc::Loc;
+use crate::maputil::{map1, map1_indexed, map2, FixedLength, FixedLengthMut};
 use crate::sets::SetOrd;
-use crate::maputil::{
-    FixedLength,
-    FixedLengthMut,
-    map1,
-    map1_indexed,
-    map2,
-};
+use crate::types::*;
+use serde::ser::{Serialize, SerializeTupleStruct, Serializer};
 
 /// A multi-dimensional cube $∏_{i=1}^N [a_i, b_i)$ with the starting and ending points
 /// along $a_i$ and $b_i$ along each dimension of type `U`.
 #[derive(Copy, Clone, Debug, Eq, PartialEq)]
-pub struct Cube<U : Num, const N : usize>(pub(super) [[U; 2]; N]);
+pub struct Cube<const N: usize, U: Num = f64>(pub(super) [[U; 2]; N]);
 
 // Need to manually implement as [F; N] serialisation is provided only for some N.
-impl<F : Num + Serialize, const N : usize> Serialize for Cube<F, N>
+impl<F: Num + Serialize, const N: usize> Serialize for Cube<N, F>
 where
     F: Serialize,
 {
@@ -50,7 +44,7 @@
     }
 }
 
-impl<A : Num, const N : usize> FixedLength<N> for Cube<A,N> {
+impl<A: Num, const N: usize> FixedLength<N> for Cube<N, A> {
     type Iter = std::array::IntoIter<[A; 2], N>;
     type Elem = [A; 2];
     #[inline]
@@ -59,7 +53,7 @@
     }
 }
 
-impl<A : Num, const N : usize> FixedLengthMut<N> for Cube<A,N> {
+impl<A: Num, const N: usize> FixedLengthMut<N> for Cube<N, A> {
     type IterMut<'a> = std::slice::IterMut<'a, [A; 2]>;
     #[inline]
     fn fl_iter_mut(&mut self) -> Self::IterMut<'_> {
@@ -67,7 +61,7 @@
     }
 }
 
-impl<'a, A : Num, const N : usize> FixedLength<N> for &'a Cube<A,N> {
+impl<'a, A: Num, const N: usize> FixedLength<N> for &'a Cube<N, A> {
     type Iter = std::slice::Iter<'a, [A; 2]>;
     type Elem = &'a [A; 2];
     #[inline]
@@ -76,15 +70,14 @@
     }
 }
 
-
 /// Iterator for [`Cube`] corners.
-pub struct CubeCornersIter<'a, U : Num, const N : usize> {
-    index : usize,
-    cube : &'a Cube<U, N>,
+pub struct CubeCornersIter<'a, U: Num, const N: usize> {
+    index: usize,
+    cube: &'a Cube<N, U>,
 }
 
-impl<'a, U : Num, const N : usize> Iterator for CubeCornersIter<'a, U, N> {
-    type Item = Loc<U, N>;
+impl<'a, U: Num, const N: usize> Iterator for CubeCornersIter<'a, U, N> {
+    type Item = Loc<N, U>;
     #[inline]
     fn next(&mut self) -> Option<Self::Item> {
         if self.index >= N {
@@ -92,28 +85,30 @@
         } else {
             let i = self.index;
             self.index += 1;
-            let arr = self.cube.map_indexed(|k, a, b| if (i>>k)&1 == 0 { a } else { b });
+            let arr = self
+                .cube
+                .map_indexed(|k, a, b| if (i >> k) & 1 == 0 { a } else { b });
             Some(arr.into())
         }
     }
 }
 
-impl<U : Num, const N : usize> Cube<U, N> {
+impl<U: Num, const N: usize> Cube<N, U> {
     /// Maps `f` over the triples $\\{(i, a\_i, b\_i)\\}\_{i=1}^N$
     /// of the cube $∏_{i=1}^N [a_i, b_i)$.
     #[inline]
-    pub fn map_indexed<T>(&self, f : impl Fn(usize, U, U) -> T) -> [T; N] {
+    pub fn map_indexed<T>(&self, f: impl Fn(usize, U, U) -> T) -> [T; N] {
         map1_indexed(self, |i, &[a, b]| f(i, a, b))
     }
 
     /// Maps `f` over the tuples $\\{(a\_i, b\_i)\\}\_{i=1}^N$
     /// of the cube $∏_{i=1}^N [a_i, b_i)$.
     #[inline]
-    pub fn map<T>(&self, f : impl Fn(U, U) -> T) -> [T; N] {
+    pub fn map<T>(&self, f: impl Fn(U, U) -> T) -> [T; N] {
         map1(self, |&[a, b]| f(a, b))
     }
 
-    /// Iterates over the start and end coordinates $\{(a_i, b_i)\}_{i=1}^N$ of the cube along 
+    /// Iterates over the start and end coordinates $\{(a_i, b_i)\}_{i=1}^N$ of the cube along
     /// each dimension.
     #[inline]
     pub fn iter_coords(&self) -> std::slice::Iter<'_, [U; 2]> {
@@ -122,27 +117,27 @@
 
     /// Returns the “start” coordinate $a_i$ of the cube $∏_{i=1}^N [a_i, b_i)$.
     #[inline]
-    pub fn start(&self, i : usize) -> U {
+    pub fn start(&self, i: usize) -> U {
         self.0[i][0]
     }
 
     /// Returns the end coordinate $a_i$ of the cube $∏_{i=1}^N [a_i, b_i)$.
     #[inline]
-    pub fn end(&self, i : usize) -> U {
+    pub fn end(&self, i: usize) -> U {
         self.0[i][1]
     }
 
     /// Returns the “start” $(a_1, … ,a_N)$ of the cube $∏_{i=1}^N [a_i, b_i)$
     /// spanned between $(a_1, … ,a_N)$ and $(b_1, … ,b_N)$.
     #[inline]
-    pub fn span_start(&self) -> Loc<U, N> {
+    pub fn span_start(&self) -> Loc<N, U> {
         Loc::new(self.map(|a, _b| a))
     }
 
     /// Returns the end $(b_1, … ,b_N)$ of the cube $∏_{i=1}^N [a_i, b_i)$
     /// spanned between $(a_1, … ,a_N)$ and $(b_1, … ,b_N)$.
     #[inline]
-    pub fn span_end(&self) -> Loc<U, N> {
+    pub fn span_end(&self) -> Loc<N, U> {
         Loc::new(self.map(|_a, b| b))
     }
 
@@ -150,19 +145,22 @@
     /// $∏_{i=1}^N [a_i, b_i)$.
     #[inline]
     pub fn iter_corners(&self) -> CubeCornersIter<'_, U, N> {
-        CubeCornersIter{ index : 0, cube : self }
+        CubeCornersIter {
+            index: 0,
+            cube: self,
+        }
     }
 
     /// Returns the width-`N`-tuple $(b_1-a_1, … ,b_N-a_N)$ of the cube $∏_{i=1}^N [a_i, b_i)$.
     #[inline]
-    pub fn width(&self) -> Loc<U, N> {
-        Loc::new(self.map(|a, b| b-a))
+    pub fn width(&self) -> Loc<N, U> {
+        Loc::new(self.map(|a, b| b - a))
     }
 
     /// Translates the cube $∏_{i=1}^N [a_i, b_i)$ by the `shift` $(s_1, … , s_N)$ to
     /// $∏_{i=1}^N [a_i+s_i, b_i+s_i)$.
     #[inline]
-    pub fn shift(&self, shift : &Loc<U, N>) -> Self {
+    pub fn shift(&self, shift: &Loc<N, U>) -> Self {
         let mut cube = self.clone();
         for i in 0..N {
             cube.0[i][0] += shift[i];
@@ -173,144 +171,158 @@
 
     /// Creates a new cube from an array.
     #[inline]
-    pub fn new(data : [[U; 2]; N]) -> Self {
+    pub fn new(data: [[U; 2]; N]) -> Self {
         Cube(data)
     }
 }
 
-impl<F : Float, const N : usize> Cube<F, N> {
+impl<F: Float, const N: usize> Cube<N, F> {
     /// Returns the centre of the cube
-    pub fn center(&self) -> Loc<F, N> {
+    pub fn center(&self) -> Loc<N, F> {
         map1(self, |&[a, b]| (a + b) / F::TWO).into()
     }
 }
 
-impl<U : Num> Cube<U, 1> {
+impl<U: Num> Cube<1, U> {
     /// Get the corners of the cube.
     ///
     /// TODO: generic implementation once const-generics can be involved in
     /// calculations.
     #[inline]
-    pub fn corners(&self) -> [Loc<U, 1>; 2] {
+    pub fn corners(&self) -> [Loc<1, U>; 2] {
         let [[a, b]] = self.0;
         [a.into(), b.into()]
     }
 }
 
-impl<U : Num> Cube<U, 2> {
+impl<U: Num> Cube<2, U> {
     /// Get the corners of the cube in counter-clockwise order.
     ///
     /// TODO: generic implementation once const-generics can be involved in
     /// calculations.
     #[inline]
-    pub fn corners(&self) -> [Loc<U, 2>; 4] {
-        let [[a1, b1], [a2, b2]]=self.0;
-        [[a1, a2].into(),
-         [b1, a2].into(),
-         [b1, b2].into(),
-         [a1, b2].into()]
+    pub fn corners(&self) -> [Loc<2, U>; 4] {
+        let [[a1, b1], [a2, b2]] = self.0;
+        [
+            [a1, a2].into(),
+            [b1, a2].into(),
+            [b1, b2].into(),
+            [a1, b2].into(),
+        ]
     }
 }
 
-impl<U : Num> Cube<U, 3> {
+impl<U: Num> Cube<3, U> {
     /// Get the corners of the cube.
     ///
     /// TODO: generic implementation once const-generics can be involved in
     /// calculations.
     #[inline]
-    pub fn corners(&self) -> [Loc<U, 3>; 8] {
-        let [[a1, b1], [a2, b2], [a3, b3]]=self.0;
-        [[a1, a2, a3].into(),
-         [b1, a2, a3].into(),
-         [b1, b2, a3].into(),
-         [a1, b2, a3].into(),
-         [a1, b2, b3].into(),
-         [b1, b2, b3].into(),
-         [b1, a2, b3].into(),
-         [a1, a2, b3].into()]
+    pub fn corners(&self) -> [Loc<3, U>; 8] {
+        let [[a1, b1], [a2, b2], [a3, b3]] = self.0;
+        [
+            [a1, a2, a3].into(),
+            [b1, a2, a3].into(),
+            [b1, b2, a3].into(),
+            [a1, b2, a3].into(),
+            [a1, b2, b3].into(),
+            [b1, b2, b3].into(),
+            [b1, a2, b3].into(),
+            [a1, a2, b3].into(),
+        ]
     }
 }
 
 // TODO: Implement Add and Sub of Loc to Cube, and Mul and Div by U : Num.
 
-impl<U : Num, const N : usize> From<[[U; 2]; N]> for Cube<U, N> {
+impl<U: Num, const N: usize> From<[[U; 2]; N]> for Cube<N, U> {
     #[inline]
-    fn from(data : [[U; 2]; N]) -> Self {
+    fn from(data: [[U; 2]; N]) -> Self {
         Cube(data)
     }
 }
 
-impl<U : Num, const N : usize> From<Cube<U, N>> for [[U; 2]; N] {
+impl<U: Num, const N: usize> From<Cube<N, U>> for [[U; 2]; N] {
     #[inline]
-    fn from(Cube(data) : Cube<U, N>) -> Self {
+    fn from(Cube(data): Cube<N, U>) -> Self {
         data
     }
 }
 
-
-impl<U, const N : usize> Cube<U, N> where U : Num + PartialOrd {
+impl<U, const N: usize> Cube<N, U>
+where
+    U: Num + PartialOrd,
+{
     /// Checks whether the cube is non-degenerate, i.e., the start coordinate
     /// of each axis is strictly less than the end coordinate.
     #[inline]
     pub fn nondegenerate(&self) -> bool {
         self.0.iter().all(|range| range[0] < range[1])
     }
-    
+
     /// Checks whether the cube intersects some `other` cube.
     /// Matching boundary points are not counted, so `U` is ideally a [`Float`].
     #[inline]
-    pub fn intersects(&self, other : &Cube<U, N>) -> bool {
-        self.iter_coords().zip(other.iter_coords()).all(|([a1, b1], [a2, b2])| {
-            a1 < b2 && a2 < b1
-        })
+    pub fn intersects(&self, other: &Cube<N, U>) -> bool {
+        self.iter_coords()
+            .zip(other.iter_coords())
+            .all(|([a1, b1], [a2, b2])| a1 < b2 && a2 < b1)
     }
 
     /// Checks whether the cube contains some `other` cube.
-    pub fn contains_set(&self, other : &Cube<U, N>) -> bool {
-        self.iter_coords().zip(other.iter_coords()).all(|([a1, b1], [a2, b2])| {
-            a1 <= a2 && b1 >= b2
-        })
+    pub fn contains_set(&self, other: &Cube<N, U>) -> bool {
+        self.iter_coords()
+            .zip(other.iter_coords())
+            .all(|([a1, b1], [a2, b2])| a1 <= a2 && b1 >= b2)
     }
 
     /// Produces the point of minimum $ℓ^p$-norm within the cube `self` for any $p$-norm.
     /// This is the point where each coordinate is closest to zero.
     #[inline]
-    pub fn minnorm_point(&self) -> Loc<U, N> {
+    pub fn minnorm_point(&self) -> Loc<N, U> {
         let z = U::ZERO;
         // As always, we assume that a ≤ b.
         self.map(|a, b| {
             debug_assert!(a <= b);
             match (a < z, z < b) {
-                (false, _)   => a,
-                (_, false)   => b,
-                (true, true) => z
+                (false, _) => a,
+                (_, false) => b,
+                (true, true) => z,
             }
-        }).into()
+        })
+        .into()
     }
 
     /// Produces the point of maximum $ℓ^p$-norm within the cube `self` for any $p$-norm.
     /// This is the point where each coordinate is furthest from zero.
     #[inline]
-    pub fn maxnorm_point(&self) -> Loc<U, N> {
+    pub fn maxnorm_point(&self) -> Loc<N, U> {
         let z = U::ZERO;
         // As always, we assume that a ≤ b.
         self.map(|a, b| {
             debug_assert!(a <= b);
             match (a < z, z < b) {
-                (false, _)   => b,
-                (_, false)   => a,
+                (false, _) => b,
+                (_, false) => a,
                 // A this stage we must have a < 0 (so U must be signed), and want to check
                 // whether |a| > |b|. We can do this without assuming U to actually implement
                 // `Neg` by comparing whether 0 > a + b.
-                (true, true) => if z > a + b { a } else { b }
+                (true, true) => {
+                    if z > a + b {
+                        a
+                    } else {
+                        b
+                    }
+                }
             }
-        }).into()
+        })
+        .into()
     }
 }
 
 macro_rules! impl_common {
     ($($t:ty)*, $min:ident, $max:ident) => { $(
-        impl<const N : usize> SetOrd for Cube<$t, N> {
+        impl<const N : usize> SetOrd for Cube<N, $t> {
             #[inline]
             fn common(&self, other : &Self) -> Self {
                 map2(self, other, |&[a1, b1], &[a2, b2]| {
@@ -338,7 +350,7 @@
 #[cfg(feature = "nightly")]
 impl_common!(f32 f64, minimum, maximum);
 
-impl<U : Num, const N : usize> std::ops::Index<usize> for Cube<U, N> {
+impl<U: Num, const N: usize> std::ops::Index<usize> for Cube<N, U> {
     type Output = [U; 2];
     #[inline]
     fn index(&self, index: usize) -> &Self::Output {
@@ -346,7 +358,7 @@
     }
 }
 
-impl<U : Num, const N : usize> std::ops::IndexMut<usize> for Cube<U, N> {
+impl<U: Num, const N: usize> std::ops::IndexMut<usize> for Cube<N, U> {
     #[inline]
     fn index_mut(&mut self, index: usize) -> &mut Self::Output {
         &mut self.0[index]