Rust, n=6 in ~12 1 second
This answer is based on @Anders Kaseorg's solution and uses some of @isaacg's improvements.
I got the answer 1,530,608,978,810 for n=7 in about 12 hours using more than 100 gigabytes of ram. (I'm not completely sure it's correct, since I don't know how to prove the conjectures this relies on).
The basic strategy is to count the total number of networks using the principle of inclusion-exclusion (splitting into groups based on the first swap).
Cargo.toml
cargo-features = ["profile-rustflags"]
[package]
name = "sorting_networks"
version = "0.1.0"
edition = "2021"
# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html
[dependencies]
hashbrown = "0.13.2"
rustc-hash = "1.1.0"
nohash-hasher = "0.2.0"
[profile.release]
rustflags = ["-C", "target-feature=+crt-static"]
main.rs
use std::{
cmp::min,
iter::repeat,
mem::{replace, swap},
time::{Duration, Instant},
};
use rustc_hash::FxHasher;
use std::hash::BuildHasherDefault;
type HashSet<V> = hashbrown::HashSet<V, BuildHasherDefault<FxHasher>>;
type HashMap<K, V> = hashbrown::HashMap<K, V, BuildHasherDefault<FxHasher>>;
macro_rules! networks {
($n:expr) => {
networks::<{ $n }, { ($n * $n - $n) / 2 }>()
};
}
type Item = u64;
fn main() {
let then = Instant::now();
networks!(1);
networks!(2);
networks!(3);
networks!(4);
networks!(5);
networks!(6);
// networks!(7);
println!("Total time: {:?}", Instant::now().duration_since(then));
}
/**
* Calculates number of networks
*/
fn networks<const N: usize, const N2: usize>() -> i64 {
let then = Instant::now();
let v: Vec<(i64, i64, Duration, Vec<u64>)> = counter(N)
.into_iter()
.map(|(sizes, multiplier)| {
let then = Instant::now();
let a = count_specific::<N, N2>(&sizes);
(a, multiplier, Instant::now().duration_since(then), sizes)
})
.collect();
let o = v.iter().map(|(a, m, _, _)| a * m).sum::<i64>() + 1;
println!(
"n={}: {} (took {:?})",
N,
o,
Instant::now().duration_since(then)
);
v.iter().for_each(|(a, m, duration, sizes)| {
println!(
" {:?}: {} = {} * {} (took {:?})",
sizes,
a * m,
a,
m,
duration
);
});
o
}
/**
* Calculate multipliers for PIE
*/
fn counter(n: usize) -> HashMap<Vec<Item>, i64> {
let mut counts = HashMap::default();
if n == 1 {
return counts;
}
let mut uf = UnionFind::new(n);
counter_(&mut uf, n, 0, 1, 1, &mut counts);
counts
}
fn counter_(
uf: &mut UnionFind,
n: usize,
i: usize,
j: usize,
m: i64,
counts: &mut HashMap<Vec<Item>, i64>,
) {
let mut j1 = j + 1;
let i1 = if j1 == n {
j1 = i + 2;
i + 1
} else {
i
};
uf.save();
uf.join(i, j);
*counts.entry(uf.get_sizes()).or_insert(0) += m;
if j1 < n {
counter_(uf, n, i1, j1, -m, counts);
}
uf.rollback();
if j1 < n {
counter_(uf, n, i1, j1, m, counts);
}
}
/**
* Counts number of networks that start with sorts that have given sizes
*/
fn count_specific<const N: usize, const N2: usize>(sizes: &Vec<Item>) -> i64 {
let mut base_state = [0; N];
let mut num_bits = 0;
if sizes.iter().filter(|&v| v & 1 == 1).count() > (N & 1) {
// no symmetric starting states
let filter: Vec<_> = sizes
.iter()
.copied()
.scan(0, |state, i| {
let offset = *state;
*state += i;
Some(((1 << i) - 1, offset))
})
.collect();
(1..1 << N)
.filter(|&v| {
// filter out states that aren't sorted in given ranges
filter.iter().all(|&(mask, offset)| {
let n = (v >> offset) & mask;
n & (n + 1) == 0
})
})
.filter(|&v: &Item| v & (v + 1) != 0) // filter out states that are completely sorted
.for_each(|v| {
// convert to state
(0..N).rev().fold(v, |acc, i| {
base_state[i] = (base_state[i] << 1) | (acc & 1);
acc >> 1
});
num_bits += 1;
});
if num_bits > Item::BITS {
panic!("Too big? {} > {} {:?}", num_bits, Item::BITS, sizes);
}
let mut states: HashSet<[Item; N]> = HashSet::default();
visit_notree::<N, N2>(&mut base_state, &mut states); // populate states
states.len() as i64
} else {
let mut state_to_idx: HashMap<[Item; N], usize> = HashMap::default();
let mut tree: Vec<[usize; N2]> = Vec::new();
let filter: Vec<_> = sizes
.iter()
.copied()
.scan(0, |state, i| {
// filter divides range so [3,2] => 0b1001001, 0b0100010
let w = i / 2;
let mask = (1 << w) - 1;
let offset0 = *state;
let mask0 = mask << offset0;
let offset1 = N as Item / 2;
let mask1 = (i & 1) << offset1;
let offset2 = N as Item - (w + offset0);
let mask2 = mask << offset2;
*state += w;
Some((
mask0,
offset0,
mask1,
offset1 - w,
mask2,
offset2 - w - (i & 1),
))
})
.collect();
(1..1 << N)
.filter(|v: &Item| {
// filter out states that aren't sorted in given ranges
filter.iter().all(|&(m0, o0, m1, o1, m2, o2)| {
let n = (v & m0) >> o0 | (v & m1) >> o1 | (v & m2) >> o2;
n & (n + 1) == 0
})
})
.filter(|&v: &Item| {
// filter out states that are completely sorted or states that have a smaller dual state
(0..=N / 2).contains(&(v.count_ones() as usize)) && (v & (v + 1) != 0)
})
.map(|v| {
// convert states to smaller states (this one is for situations like 010011 vs 001101)
min(
v,
(((1 << N) - 1) - v).reverse_bits() >> (Item::BITS as usize - N),
)
})
.collect::<HashSet<Item>>() //remove duplicates
.into_iter()
.for_each(|v| {
// convert to state
(0..N).rev().fold(v, |acc, i| {
base_state[i] = (base_state[i] << 1) | (acc & 1);
acc >> 1
});
num_bits += 1;
});
if num_bits > Item::BITS {
panic!("Too big? {} > {} {:?}", num_bits, Item::BITS, sizes);
} // currently, the largest this gets (for sizes=[2] and N=7) is 44
let r = visit(&mut base_state, &mut tree, &mut state_to_idx); // build tree ("tree")
double::<N, N2>(tree, r) as i64 // doubles (more like squares) tree
}
}
/**
* Finds all states that are descendants of the input state and maps them to numbers in topological order.
* "tree" is an adjacency list on the numbers
*/
fn visit<const N: usize, const N2: usize>(
state: &mut [Item; N],
tree: &mut Vec<[usize; N2]>,
state_to_idx: &mut HashMap<[Item; N], usize>,
) -> usize {
if let Some(&idx) = state_to_idx.get(state) {
idx
} else {
let mut l: Vec<Option<usize>> = vec![];
for a in 0..N - 1 {
let p = state[a];
for b in a + 1..N {
let q = state[b];
if p & !q == 0 {
// self loop
l.push(None);
} else {
state[a] = p & q;
state[b] = p | q;
l.push(Some(visit::<N, N2>(state, tree, state_to_idx)));
state[a] = p;
state[b] = q;
}
}
}
let v = state_to_idx.len();
state_to_idx.insert(*state, v);
tree.push(
l.iter()
.copied()
.map(|x| x.unwrap_or(v))
.collect::<Vec<usize>>()
.try_into()
.unwrap(),
);
v
}
}
/**
* Finds all descendant states from input_state
*/
fn visit_notree<const N: usize, const N2: usize>(
state: &mut [Item; N],
states: &mut HashSet<[Item; N]>,
) {
if !states.contains(state) {
for a in 0..N - 1 {
let p = state[a];
for b in a + 1..N {
let q = state[b];
if p & !q != 0 {
state[a] = p & q;
state[b] = p | q;
visit_notree::<N, N2>(state, states);
state[a] = p;
state[b] = q;
}
}
}
states.insert(*state);
}
}
/**
* Combines a graph with the reverse of the graph
*/
fn double<const N: usize, const N2: usize>(tree: Vec<[usize; N2]>, r: usize) -> usize {
// does bfs on the graph, because visit returns in topological order, we can discard states with smaller rank, which saves memory
let choices: Vec<(usize, usize)> = (0..N - 1).flat_map(|a| repeat(a).zip(a + 1..N)).collect();
let choices_idx: HashMap<(usize, usize), usize> = choices
.iter()
.copied()
.enumerate()
.map(|(a, b)| (b, a))
.collect();
let rev: Vec<usize> = choices
.iter()
.copied()
.map(|(a, b)| *choices_idx.get(&(N - 1 - b, N - 1 - a)).unwrap())
.collect();
let tree_rev: Vec<[usize; N2]> = tree
.iter()
.map(|t| {
rev.iter()
.map(|&v| t[v])
.collect::<Vec<usize>>()
.try_into()
.unwrap()
})
.collect();
let tl = tree.len();
let mut queues: Vec<HashSet<usize>> = (0..2 * tl - 1).map(|_| HashSet::default()).collect();
queues[2 * r].insert(r);
let mut center = 0;
let mut num_seen = 0;
let mut rank = 2 * tl - 1;
while !queues.is_empty() {
rank -= 1;
queues.pop().unwrap().drain().for_each(|a| {
let b = rank - a;
if a == b {
center += 1;
}
num_seen += 1;
for (&a1, &b1) in tree[a].iter().zip(tree_rev[b].iter()) {
if a1 > b1 {
if a != b && (a != b1 || b != a1) {
queues[b1 + a1].insert(b1);
}
} else if a != a1 || b != b1 {
queues[a1 + b1].insert(a1);
}
}
});
}
2 * num_seen - center
}
// union-find with backtracking whatever
struct BacktrackArray<T> {
data: Vec<T>,
history: Vec<(usize, T)>,
checkpoints: Vec<usize>,
}
impl<T: Copy> BacktrackArray<T> {
pub fn new(data: Vec<T>) -> BacktrackArray<T> {
BacktrackArray {
data,
history: vec![],
checkpoints: vec![],
}
}
pub fn rollback(&mut self) {
self.history
.drain(self.checkpoints.pop().unwrap_or(0)..)
.rev()
.for_each(|(idx, v)| self.data[idx] = v);
}
pub fn save(&mut self) {
self.checkpoints.push(self.history.len());
}
pub fn set(&mut self, idx: usize, value: T) {
self.history
.push((idx, replace(&mut self.data[idx], value)));
}
pub fn get(&self, idx: usize) -> T {
self.data[idx]
}
}
impl<V: Copy> FromIterator<V> for BacktrackArray<V> {
fn from_iter<T: IntoIterator<Item = V>>(iter: T) -> Self {
BacktrackArray::new(Vec::from_iter(iter))
}
}
impl<V: Copy> From<Vec<V>> for BacktrackArray<V> {
fn from(val: Vec<V>) -> Self {
BacktrackArray::new(val)
}
}
struct UnionFind {
num_unions: BacktrackArray<usize>,
parent: BacktrackArray<usize>,
size: BacktrackArray<Item>,
representatives: BacktrackArray<usize>,
rep_ptrs: BacktrackArray<usize>,
}
impl UnionFind {
pub fn new(n: usize) -> UnionFind {
UnionFind {
num_unions: vec![n].into(),
parent: (0..n).collect(),
size: repeat(1).take(n).collect(),
representatives: (0..n).collect(),
rep_ptrs: (0..n).collect(),
}
}
pub fn save(&mut self) {
self.num_unions.save();
self.parent.save();
self.size.save();
self.representatives.save();
self.rep_ptrs.save();
}
pub fn rollback(&mut self) {
self.num_unions.rollback();
self.parent.rollback();
self.size.rollback();
self.representatives.rollback();
self.rep_ptrs.rollback();
}
pub fn get_sizes(&self) -> Vec<Item> {
let n = self.num_unions.get(0);
let mut v: Vec<Item> = self
.representatives
.data
.iter()
.take(n)
.map(|&i| self.size.get(i))
.filter(|&v| v > 1)
.collect();
v.sort_unstable();
v
}
pub fn join(&mut self, i: usize, j: usize) {
let mut i = self.find(i);
let mut j = self.find(j);
if i != j {
let mut s_i = self.size.get(i);
let mut s_j = self.size.get(j);
if s_i < s_j {
swap(&mut s_i, &mut s_j);
swap(&mut i, &mut j);
}
self.parent.set(j, i);
self.size.set(i, s_i + s_j);
let n = self.num_unions.get(0) - 1;
self.num_unions.set(0, n);
let r = self.representatives.get(n);
let r_idx = self.rep_ptrs.get(j);
self.representatives.set(r_idx, r);
self.rep_ptrs.set(r, r_idx);
}
}
pub fn find(&mut self, i: usize) -> usize {
let i1 = self.parent.get(i);
if i1 == i {
i
} else {
let i2 = self.find(i1);
self.parent.set(i, i2);
i2
}
}
}
Output:
n=1: 1
n=2: 2
[2]: 1 = 1 * 1 (took 7.73µs)
n=3: 11
[2]: 12 = 4 * 3 (took 9.617µs)
[3]: -2 = 1 * -2 (took 783ns)
n=4: 261
[2]: 366 = 61 * 6 (took 17.105µs)
[2, 2]: -48 = 16 * -3 (took 3.853µs)
[3]: -64 = 8 * -8 (took 2.727µs)
[4]: 6 = 1 * 6 (took 1.088µs)
n=5: 43337
[2]: 62480 = 6248 * 10 (took 440.742µs)
[2, 2]: -14100 = 940 * -15 (took 76.582µs)
[5]: -24 = 1 * -24 (took 1.818µs)
[2, 3]: 1280 = 64 * 20 (took 9.063µs)
[3]: -6780 = 339 * -20 (took 35.228µs)
[4]: 480 = 16 * 30 (took 3.91µs)
n=6: 72462128
[2, 3]: 1741440 = 14512 * 120 (took 3.553571ms)
[4]: 169650 = 1885 * 90 (took 164.134µs)
[2, 2, 2]: 831870 = 55458 * 15 (took 3.841656ms)
[2]: 103912905 = 6927527 * 15 (took 553.271018ms)
[2, 2]: -28154610 = 625658 * -45 (took 43.410737ms)
[3, 3]: -20480 = 512 * -40 (took 67.739µs)
[5]: -4608 = 32 * -144 (took 4.743µs)
[6]: 120 = 1 * 120 (took 2.899µs)
[3]: -5991120 = 149778 * -40 (took 33.794472ms)
[2, 4]: -23040 = 256 * -90 (took 49.235µs)
n=7: 1530608978810
[2, 2]: -701933448825 = 6685080465 * -105 (took 1777.016031993s)
[7]: -720 = 1 * -720 (took 22.07µs)
[2, 2, 3]: -684133800 = 3257780 * -210 (took 424.2515ms)
[2, 3]: 26000250360 = 61905358 * 420 (took 7.892396883s)
[4]: 754063170 = 3590777 * 210 (took 282.065867ms)
[2, 2, 2]: 39011470995 = 371537819 * 105 (took 67.801355981s)
[2]: 2244447515625 = 106878453125 * 21 (took 40712.141796623s)
[3, 4]: 1720320 = 4096 * 420 (took 431.995µs)
[2, 5]: 516096 = 1024 * 504 (took 93.822µs)
[5]: -5283432 = 10483 * -504 (took 809.017µs)
[6]: 53760 = 64 * 840 (took 13.258µs)
[3]: -76668548810 = 1095264983 * -70 (took 215.344123884s)
[3, 3]: -174011040 = 621468 * -280 (took 176.875932ms)
[2, 4]: -141184890 = 224103 * -630 (took 93.090671ms)
Total time: 42782.049221178s
The correctness of this solution relies on two conjectures.
Before that, notation.
There are \$N\$ wires, \$1\$ through \$N\$ (the set of them is denoted by \$[N]\$).
\$U\$ is the set of all possible networks over the \$N\$ wires. \$2^U\$ is the powerset of \$U\$ (that is, the set of all subsets of \$U\$).
We'll let \$G=\ (U, E)\$ be a directed graph with vertices \$U\$ and an edge from \$a\$ to \$b\$ iff \$b\$ can be expressed by appending a single vertical wire to the right of \$a\$. We also denote by \$G_x\$ (for network \$x\$) the set of networks reachable from \$x\$.
Define \$S(A):\ 2^{[N]}\rightarrow U\$ to mean the network that sorts the wires in \$A\$.
This is sort of an abuse of notation, but for \$A\subseteq[N]\$, we will write \$G_{S(A)}\$ as \$G_A\$
Define \$F(\{A_1, ..., A_k\}): 2^{\left(2^{[N]}\right)}\rightarrow 2^U\$ (where \$A_i \subseteq [N]\$ and the \$A_i\$ are pairwise disjoint) to mean \$G_{A_1}\cap G_{A_2}\cap\ ...\cap\ G_{A_k}\$.
The first conjecture is that \$\|F(\{A_1,...,A_k\})\|=\|F(\{B_1,...,B_k\})\|\$ if \$\|A_i\|=\|B_i\|\$ for \$1\le i\le k\$.
The second conjecture is that given \$X_1, X_2, ..., X_m \subseteq [N]\$ (possibly overlapping), \$G_{X_1}\cap G_{X_2}\cap ...\cap\ G_{X_m}=F({A_1,A_2,...,A_k})\$ for some pairwise disjoint \$A_i\$ with \$k\le m\$.
To be more specific, if we make a graph with vertices in \$[N]\$ and an edge between vertices \$a\$ and \$b\$ if there is some \$i\$ for which \$a,b \in X_i\$, then the \$A_j\$ are the connected components of that graph.
