[WIP] Parallel linclust

[bbuschkaemper]

Shared-filesystem parallel linclust

Work in progress. Runs linclust across many Slurm workers on a shared filesystem, with no MPI and
no node-to-node communication, to allow clustering 1e11-1e12 sequences.

Why

We assume availability of up to 2TB memory nodes and large amounts of shared filesystem storage ("scratch").
Several structures in linclust are sized by key space or by the whole database, so they cannot
exist at 1e12 on a 2 TB node:

structure	at 1e12
seqkey_to_len, countTable, repSequence (kmermatcher.cpp)	2-4 TB each
resident DBReader::Index[]	~24 TB
assignedCluster[dbSize] (Align2clust.cpp:445)	8 TB
std::list<size_t>[N] (mergeclusters.cpp:28)	24 TB of empty headers

There is also a time cost. When the k-mer array does not fit memory, kmermatcher splits and
re-extracts every k-mer per split, roughly 68 times over at 1e12.

Approach

Partition k-mer space rather than sequence space. The partition of a k-mer is the low bits of the
hashUInt64 score kmermatcher already computes, so every occurrence of a k-mer lands in the same
partition and a partition can be grouped on its own. The division is lossless. Sequence-space
sharding is not.

Keys are dense and length-ranked: createdbparallel writes sequences longest first, so a key is
its global length rank. That removes SORT_BY_LENGTH and its side arrays, lets an entry be
addressed by key with no resident index, and makes stock's longest-first greedy the same thing as
ascending key order. Clustering becomes one left-to-right sweep needing 2 bits per key instead of
8 bytes, so 25 GB at 1e11 rather than 800 GB.

Coordination is files only. Every worker of a stage runs the same command line and takes its
identity from a fetch_add on a counter file, so a stage maps onto a Slurm array job and workers
can join late, die, or restart. Items are claimed under a lease and re-claimed if a worker dies.
Locking is fcntl whole-file locks, not flock, which is node-local on GPFS and Lustre.

Nine commands, driven by data/workflow/linclustparallel.sh:

createdbparallel        FASTA -> dense, length-ranked sequence DB
kmermatcherparallel     extract k-mers into P partition buckets, one wave at a time
kmerreduceparallel      per partition: group, emit candidate edges bucketed by representative key
alignparallel           per bucket: merge duplicate copies of a pair, align once
greedycluster           single node, one key-ordered sweep
createrepdb / translatecluster / mergeclusterparallel / translatekeys

Two decisions that were measured, not assumed

Alignment is keyed by representative range, not by k-mer partition. We built the fused version
first. A k-mer partition's pairs are spread over the whole key space, so one partition needed
1.61 GB of sequences and read 83.9 GB, a 52x amplification. Bucketing edges by representative key
gives 1.15x, removes cross-partition duplicates before aligning rather than paying for them, and
reproduces stock's global per-(pair, diagonal) accumulation exactly.

Bucket records are packed. K-mer records went from 24 to ~14 bytes, candidate edges from 17 to
~7. --raw-records writes the old fixed-width form, so the two can be run against each other;
they produce identical output. Peak scratch is 0.67x what the first version of this branch used.

Results

Verified identical results on a 1M subsample of MGnify sequences, current master branch ("stock") and this branch on the same database:

0 of 1,000,000 sequences in a differing cluster
seq -> rep identical for 1,000,000 of 1,000,000

Output is byte-identical across worker counts, wave counts, thread counts, record encodings, reduce slice counts, 32- and 64-bit builds, and after kill-and-resume.

Comparing against stock on a differently keyed database is not meaningful: stock against itself
on input-order versus length-ranked keys moves 3.83% of sequences.

Speed on one machine, 128 cores total, FASTA to TSV:

	stock	this branch, 8 workers
1M	21 s	50 s
10M	137 s	121 s
100M	1334 s	1078 s

Peak memory at 100M is 61.9 GB against stock's 91.6 GB. On a 64-bit-id build, which a real run
needs, the gap widens: 5.7 GB against 11.2 GB at 10M, because stock's key-space-sized structures
scale with key width and these do not.

Splitting the same cores further has diminishing returns: 1 to 4 workers is worth 2.7x, 4 to 8
only 1.08x.

Not finished

• Multi-node, real slurm was verified on smaller subsamples (1-10M) only, at least a 1B stock vs. multi-node parallel run should be done.
• At 1e12 the map has ~1e6 work items, each costing two fsyncs through one global lock. This needs batched claims and sharded queue files.
• No checks up front that linclust's pass-2 will fit the scratch budget, so a long run can still fail late.

Limitations

• --cov-mode 1 or 2 only. Symmetric coverage modes make linclust select SET_COVER plus the count-table rounds; neither is implemented, so the command refuses them.
• Protein only. alignparallel rejects nucleotides.
• Output is representative<TAB>member in accessions, not a cluster DB. A per-key index is state no single node can hold at this scale.
• Sequences are capped at 65,535 residues and rejected above it.
• Stock is touched in 6 files (+150/−77): two defaulted-NULL parameters and guarded branches in
  kmermatcher.{cpp,h}, and parsePrecisionLib de-duplicated into Matcher.cpp. Both were
  checked behaviour-preserving against a reference binary built from the base commit.
  Parameters.cpp also relaxes checkIfDatabaseIsValid so a mkdir race is tolerated when the
  directory already exists, which fixes a real Slurm-array race.