diff --git a/CLAUDE.md b/CLAUDE.md index ee5c865..90e036f 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -42,9 +42,13 @@ Uses standard double-dash flags (`--contig`, `--abund`, etc.). Single-dash compa - `maxbin-rs filter` — filter contigs by minimum length - `maxbin-rs seeds` — generate seed file from HMMER marker gene hits - `maxbin-rs em` — run the Rust EM algorithm -- `maxbin-rs cpp-em` — run the original C++ EM via FFI (equivalence testing only) - `maxbin-rs sam-to-abund` — compute abundance from a SAM file +The original C++ EM (used only for equivalence testing) now lives in a +separate crate: run it via `maxbin-rs-cpp-em` from `crates/maxbin-rs-equivalence` +(see below). It is no longer a `maxbin-rs` subcommand, so the core crate carries +no C/C++. + ### Key source files (in `crates/maxbin-rs/src/`) - `pipeline.rs` — orchestration for each subcommand @@ -53,10 +57,30 @@ Uses standard double-dash flags (`--contig`, `--abund`, etc.). Single-dash compa - `kmer_map.rs` — tetranucleotide frequency computation - `abundance.rs` — abundance file parsing - `fasta.rs` — FASTA reading/writing -- `original_ffi.rs` — C++ FFI bridge for equivalence testing against the original EM - `external.rs` — shelling out to HMMER, Bowtie2, gene caller - `profiler.rs` — timing instrumentation +The core crate has no `build.rs` and no C/C++ dependency, so it can eventually +compile to wasm. + +### Equivalence crate (`crates/maxbin-rs-equivalence/`) + +Holds everything that links the original MaxBin2 C++ for FFI-based equivalence +testing, kept out of the core crate: + +- `build.rs` + `vendor/`: extract the upstream MaxBin2 C++ from the tarball + (`MAXBIN2_SRC_TARBALL`), apply `vendor/maxbin2-cpp-ffi.patch`, and compile it + as a static lib via `cc`. +- `src/original_ffi.rs`: Rust FFI bindings to that C++ (re-exported from the + crate's `lib.rs`). +- `src/bin/maxbin-rs-cpp-em.rs`: the `maxbin-rs-cpp-em` binary that runs the + original C++ EM via FFI (formerly the `cpp-em` subcommand). +- `tests/`: the FFI equivalence tests and the FFI-comparing proptests. They + exercise the core crate's public API and compare it against the C++. The + `divergent-em` fixture lives here too. + +The crate depends on `maxbin-rs` as a path dependency. + ### Parallelism Both the inner abundance-file loop and the outer contig loop in the EM E-step are parallelized with Rayon, controlled by `--thread`. @@ -74,7 +98,7 @@ Outputs: one FASTA per bin + summary statistics + `.log` file. ## Testing Strategy -- **Equivalence tests** compare Rust output against original MaxBin2 on the same input, including reproducing known bugs. The C++ EM is available via FFI (`cpp-em` subcommand) for component-level comparison. +- **Equivalence tests** compare Rust output against original MaxBin2 on the same input, including reproducing known bugs. The C++ EM is available via FFI in the `maxbin-rs-equivalence` crate (the `maxbin-rs-cpp-em` binary) for component-level comparison. - **Property tests** (proptest) cover parsing, distance computation, EM numerics, and sorting. - **Benchmarks** via `cargo nextest run --cargo-profile bench bench_components -- --ignored`. - When an equivalence test reveals a bug in the original, **document it** (in TODO.md and/or inline) but still reproduce the buggy behavior in the v0.1.x code path. diff --git a/Cargo.lock b/Cargo.lock index 4a64eb4..ff02fb3 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -430,11 +430,9 @@ checksum = "5e5032e24019045c762d3c0f28f5b6b8bbf38563a65908389bf7978758920897" name = "maxbin-rs" version = "0.4.0" dependencies = [ - "cc", "clap 4.6.1", "flate2", "frag_gene_scan_rs", - "proptest", "rayon", "seq_io", "serde", @@ -442,6 +440,17 @@ dependencies = [ "tempfile", ] +[[package]] +name = "maxbin-rs-equivalence" +version = "0.4.0" +dependencies = [ + "cc", + "clap 4.6.1", + "flate2", + "maxbin-rs", + "proptest", +] + [[package]] name = "memchr" version = "2.8.0" diff --git a/Cargo.toml b/Cargo.toml index 71838e4..f29d547 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -1,6 +1,6 @@ [workspace] resolver = "3" -members = ["crates/maxbin-rs"] +members = ["crates/maxbin-rs", "crates/maxbin-rs-equivalence"] # Profiles are only honored at the workspace root, so the benchmark profile # that used to live in the crate manifest moves here. diff --git a/crates/maxbin-rs-equivalence/Cargo.toml b/crates/maxbin-rs-equivalence/Cargo.toml new file mode 100644 index 0000000..0afbabd --- /dev/null +++ b/crates/maxbin-rs-equivalence/Cargo.toml @@ -0,0 +1,16 @@ +[package] +name = "maxbin-rs-equivalence" +version = "0.4.0" +edition = "2024" +publish = false + +[dependencies] +clap = { version = "4.6.1", features = ["derive"] } +maxbin-rs = { path = "../maxbin-rs" } + +[build-dependencies] +cc = "1.2.60" + +[dev-dependencies] +flate2 = "1.1.9" +proptest = "1.11.0" diff --git a/crates/maxbin-rs/build.rs b/crates/maxbin-rs-equivalence/build.rs similarity index 100% rename from crates/maxbin-rs/build.rs rename to crates/maxbin-rs-equivalence/build.rs diff --git a/crates/maxbin-rs-equivalence/src/bin/maxbin-rs-cpp-em.rs b/crates/maxbin-rs-equivalence/src/bin/maxbin-rs-cpp-em.rs new file mode 100644 index 0000000..e8e5596 --- /dev/null +++ b/crates/maxbin-rs-equivalence/src/bin/maxbin-rs-cpp-em.rs @@ -0,0 +1,48 @@ +//! `maxbin-rs-cpp-em`: run the original MaxBin2 C++ EM via FFI. +//! +//! This binary exists only for equivalence testing and benchmarking against +//! the original implementation. The core `maxbin-rs` crate no longer carries +//! the C++ FFI, so this lives in the equivalence crate. +//! +//! The C++ EManager handles its own data loading, EM, classification, and +//! file writing internally, so the whole thing is timed as one block. +use std::path::PathBuf; + +use clap::Parser; +use maxbin_rs_equivalence::original_ffi::OriginalEManager; + +#[derive(Parser)] +#[command( + name = "maxbin-rs-cpp-em", + about = "Run the original MaxBin2 C++ EM via FFI (equivalence testing)." +)] +struct Args { + #[arg(long)] + contig: PathBuf, + #[arg(long)] + abund: PathBuf, + #[arg(long)] + seed: PathBuf, + #[arg(long)] + out: String, + #[arg(long, default_value_t = 1)] + thread: usize, +} + +fn main() { + let args = Args::parse(); + + eprintln!("Running original C++ EM via FFI..."); + eprintln!(" contigs: {}", args.contig.display()); + eprintln!(" abund: {}", args.abund.display()); + eprintln!(" seed: {}", args.seed.display()); + eprintln!(" threads: {}", args.thread); + + let t = std::time::Instant::now(); + let em = OriginalEManager::new(&args.contig, &args.abund, &args.out); + em.set_thread_num(args.thread as i32); + let result = em.run(&args.seed); + let secs = t.elapsed().as_secs_f64(); + + eprintln!("C++ EM completed in {secs:.1}s (returned: {result})"); +} diff --git a/crates/maxbin-rs-equivalence/src/lib.rs b/crates/maxbin-rs-equivalence/src/lib.rs new file mode 100644 index 0000000..1668a4d --- /dev/null +++ b/crates/maxbin-rs-equivalence/src/lib.rs @@ -0,0 +1,14 @@ +//! Equivalence-testing harness for maxbin-rs. +//! +//! This crate links the original MaxBin2 C++ code as a static library and +//! exposes Rust FFI bindings to it (`original_ffi`). It exists so the core +//! `maxbin-rs` crate can stay free of any C/C++ toolchain dependency (and +//! eventually compile to wasm), while the bit-for-bit equivalence tests +//! against the original implementation keep running here. +//! +//! The C++ sources are extracted from the upstream MaxBin2 tarball and +//! patched for FFI at build time (see `build.rs` and `vendor/`). The +//! `MAXBIN2_SRC_TARBALL` environment variable points at that tarball; the +//! Nix devshell sets it automatically. + +pub mod original_ffi; diff --git a/crates/maxbin-rs/src/original_ffi.rs b/crates/maxbin-rs-equivalence/src/original_ffi.rs similarity index 99% rename from crates/maxbin-rs/src/original_ffi.rs rename to crates/maxbin-rs-equivalence/src/original_ffi.rs index e52f22e..bd28590 100644 --- a/crates/maxbin-rs/src/original_ffi.rs +++ b/crates/maxbin-rs-equivalence/src/original_ffi.rs @@ -1,5 +1,5 @@ /// FFI bindings to the original MaxBin2 C++ code. -/// Used only for equivalence testing — not compiled into release builds. +/// Used only for equivalence testing, inside the maxbin-rs-equivalence crate. use std::ffi::{CStr, CString}; use std::os::raw::c_char; diff --git a/crates/maxbin-rs-equivalence/tests/abundance_ffi.rs b/crates/maxbin-rs-equivalence/tests/abundance_ffi.rs new file mode 100644 index 0000000..106bc03 --- /dev/null +++ b/crates/maxbin-rs-equivalence/tests/abundance_ffi.rs @@ -0,0 +1,43 @@ +//! FFI equivalence test for the abundance-file parser. +//! +//! Compares `maxbin_rs::abundance` against the original MaxBin2 C++ +//! AbundanceLoader. Previously an in-source test in `abundance.rs`. +use maxbin_rs::abundance::parse_file; +use maxbin_rs_equivalence::original_ffi::OriginalAbundanceLoader; + +#[test] +fn equivalence_with_original() { + // Writes a temp abundance file, parses with both Rust and C++ FFI, + // compares field by field. + use std::io::Write; + + let test_data = "contig_1\t10.5\ncontig_2\t0.00001\ncontig_3\t100.0\ncontig_4 0.5\n"; + + let tmp_path = std::env::temp_dir().join("maxbin_rs_test_abund.txt"); + { + let mut f = std::fs::File::create(&tmp_path).unwrap(); + f.write_all(test_data.as_bytes()).unwrap(); + } + + // Rust + let rust_records = parse_file(&tmp_path).unwrap(); + + // Original C++ + let original = OriginalAbundanceLoader::new(&tmp_path); + + assert_eq!(rust_records.len(), original.num_records() as usize); + assert!(original.is_parse_success()); + + for (i, rec) in rust_records.iter().enumerate() { + let orig_abund = original.abundance_by_index(i as i32); + assert!( + (rec.abundance - orig_abund).abs() < 1e-10, + "abundance mismatch at record {i} (header: {}): rust={} original={}", + rec.header, + rec.abundance, + orig_abund, + ); + } + + std::fs::remove_file(&tmp_path).ok(); +} diff --git a/crates/maxbin-rs/tests/bench_components.rs b/crates/maxbin-rs-equivalence/tests/bench_components.rs similarity index 93% rename from crates/maxbin-rs/tests/bench_components.rs rename to crates/maxbin-rs-equivalence/tests/bench_components.rs index 6f41b4a..c62b798 100644 --- a/crates/maxbin-rs/tests/bench_components.rs +++ b/crates/maxbin-rs-equivalence/tests/bench_components.rs @@ -67,7 +67,8 @@ fn bench_normal_distribution() -> (f64, f64) { let inputs: Vec = (0..1000).map(|i| i as f64 * 0.001).collect(); let cpp_count = count_ops(BENCH_DURATION, || { - let nd = maxbin_rs::original_ffi::OriginalNormalDistribution::new(mean, std_dev); + let nd = + maxbin_rs_equivalence::original_ffi::OriginalNormalDistribution::new(mean, std_dev); for &x in &inputs { std::hint::black_box(nd.prob(x)); } @@ -92,7 +93,7 @@ fn bench_quicksort() -> (f64, f64) { let cpp_count = count_ops(BENCH_DURATION, || { let mut arr = base.clone(); let mut idx: Vec = (0..arr.len() as i32).collect(); - maxbin_rs::original_ffi::original_quicksort(&mut arr, Some(&mut idx)); + maxbin_rs_equivalence::original_ffi::original_quicksort(&mut arr, Some(&mut idx)); std::hint::black_box(&arr); }); @@ -125,7 +126,7 @@ fn bench_kmermap_lookup() -> (f64, f64) { }; let cpp_count = count_ops(BENCH_DURATION, || { - let km = maxbin_rs::original_ffi::OriginalKmerMap::new(4, true); + let km = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(4, true); for kmer in &kmers { std::hint::black_box(km.get_mapping(kmer)); } @@ -145,10 +146,10 @@ fn bench_kmermap_lookup() -> (f64, f64) { fn bench_profiler(seqs: &[Vec]) -> (f64, f64) { let cpp_count = count_ops(BENCH_DURATION, || { - let km = maxbin_rs::original_ffi::OriginalKmerMap::new(4, true); + let km = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(4, true); for seq in seqs { let s = std::str::from_utf8(seq).unwrap(); - let prof = maxbin_rs::original_ffi::OriginalProfiler::new(4, s, &km); + let prof = maxbin_rs_equivalence::original_ffi::OriginalProfiler::new(4, s, &km); std::hint::black_box(prof.get_profile(136)); } }); @@ -168,7 +169,7 @@ fn bench_profiler(seqs: &[Vec]) -> (f64, f64) { fn bench_euc_dist(profiles: &[(Vec, Vec)]) -> (f64, f64) { let cpp_count = count_ops(BENCH_DURATION, || { - let cpp = maxbin_rs::original_ffi::OriginalEucDist::new(4); + let cpp = maxbin_rs_equivalence::original_ffi::OriginalEucDist::new(4); for (p1, p2) in profiles { std::hint::black_box(cpp.get_dist_profile(p1, p2)); } @@ -187,7 +188,7 @@ fn bench_euc_dist(profiles: &[(Vec, Vec)]) -> (f64, f64) { fn bench_spearman_dist(profiles: &[(Vec, Vec)]) -> (f64, f64) { let cpp_count = count_ops(BENCH_DURATION, || { - let cpp = maxbin_rs::original_ffi::OriginalSpearmanDist::new(4); + let cpp = maxbin_rs_equivalence::original_ffi::OriginalSpearmanDist::new(4); for (p1, p2) in profiles { std::hint::black_box(cpp.get_dist_profile(p1, p2)); } @@ -207,7 +208,7 @@ fn bench_spearman_dist(profiles: &[(Vec, Vec)]) -> (f64, f64) { fn bench_fasta_parse(path: &std::path::Path, raw_bytes: &[u8]) -> (f64, f64) { let cpp_count = count_ops(BENCH_DURATION, || { - let reader = maxbin_rs::original_ffi::OriginalFastaReader::new(path); + let reader = maxbin_rs_equivalence::original_ffi::OriginalFastaReader::new(path); std::hint::black_box(reader.num_records()); }); @@ -223,7 +224,7 @@ fn bench_fasta_parse(path: &std::path::Path, raw_bytes: &[u8]) -> (f64, f64) { fn bench_abundance_parse(path: &std::path::Path, raw_bytes: &[u8]) -> (f64, f64) { let cpp_count = count_ops(BENCH_DURATION, || { - let loader = maxbin_rs::original_ffi::OriginalAbundanceLoader::new(path); + let loader = maxbin_rs_equivalence::original_ffi::OriginalAbundanceLoader::new(path); std::hint::black_box(loader.num_records()); }); diff --git a/crates/maxbin-rs-equivalence/tests/distance_ffi.rs b/crates/maxbin-rs-equivalence/tests/distance_ffi.rs new file mode 100644 index 0000000..620a06a --- /dev/null +++ b/crates/maxbin-rs-equivalence/tests/distance_ffi.rs @@ -0,0 +1,144 @@ +//! Component-level FFI equivalence tests for distance metrics. +//! +//! These compare the Rust distance functions in `maxbin_rs::distance` against +//! the original MaxBin2 C++ EucDist / SpearmanDist via FFI. They previously +//! lived as in-source `#[cfg(test)]` tests inside `distance.rs`; they moved +//! here when the C++ FFI was extracted from the core crate. +use maxbin_rs::distance::{ + DistanceContext, euc_dist_profiles, euc_dist_seq, spearman_dist_profiles, spearman_dist_seq, +}; +use maxbin_rs::kmer_map::KmerMap; +use maxbin_rs::profiler::Profiler; +use maxbin_rs_equivalence::original_ffi::{OriginalEucDist, OriginalSpearmanDist}; + +#[test] +fn ffi_euc_dist_seq_equivalence() { + let test_pairs = [ + ("ACGTACGTACGTACGT", "ACGTACGTACGTACGT"), + ("AAAAAAAAAAAAAAAA", "CCCCCCCCCCCCCCCC"), + ("ATCGATCGATCGATCG", "TAGCTAGCTAGCTAGC"), + ("ACGTACGTACGTACGTACGTACGT", "TTTTAAAACCCCGGGG"), + ("AACCGGTTAACCGGTT", "GGTTAACCGGTTAACC"), + ]; + + let ctx = DistanceContext::new(4); + let cpp = OriginalEucDist::new(4); + + for (seq1, seq2) in &test_pairs { + let rust_d = euc_dist_seq(&ctx, seq1.as_bytes(), seq2.as_bytes()); + let cpp_d = cpp.get_dist_seq(seq1, seq2); + + assert_eq!( + rust_d.to_bits(), + cpp_d.to_bits(), + "EucDist seq not bit-identical for '{seq1}' vs '{seq2}': rust={rust_d:e} cpp={cpp_d:e}" + ); + } +} + +#[test] +fn ffi_euc_dist_profile_equivalence() { + // Build profiles from known sequences, then compare profile-based distance + let cpp = OriginalEucDist::new(4); + + let seqs = ["ACGTACGTACGTACGT", "AAAAAAAAAAAAAAAA", "ATCGATCGATCGATCG"]; + + let rust_kmap = KmerMap::new(4, true); + let profiles: Vec> = seqs + .iter() + .map(|s| { + Profiler::new(4, s.as_bytes(), &rust_kmap) + .get_profile() + .to_vec() + }) + .collect(); + + for i in 0..profiles.len() { + for j in (i + 1)..profiles.len() { + let rust_d = euc_dist_profiles(&profiles[i], &profiles[j]); + let cpp_d = cpp.get_dist_profile(&profiles[i], &profiles[j]); + + assert_eq!( + rust_d.to_bits(), + cpp_d.to_bits(), + "EucDist profile not bit-identical: rust={rust_d:e} cpp={cpp_d:e}" + ); + } + } +} + +#[test] +fn ffi_spearman_dist_seq_equivalence() { + let test_pairs = [ + ("ACGTACGTACGTACGT", "ACGTACGTACGTACGT"), + ("AAAAAAAAAAAAAAAA", "CCCCCCCCCCCCCCCC"), + ("ATCGATCGATCGATCG", "TAGCTAGCTAGCTAGC"), + ("ACGTACGTACGTACGTACGTACGT", "TTTTAAAACCCCGGGG"), + ("AACCGGTTAACCGGTT", "GGTTAACCGGTTAACC"), + ]; + + let ctx = DistanceContext::new(4); + let cpp = OriginalSpearmanDist::new(4); + + for (seq1, seq2) in &test_pairs { + let rust_d = spearman_dist_seq(&ctx, seq1.as_bytes(), seq2.as_bytes()); + let cpp_d = cpp.get_dist_seq(seq1, seq2); + + let diff = (rust_d - cpp_d).abs(); + assert!( + diff < 1e-10, + "SpearmanDist seq mismatch for '{seq1}' vs '{seq2}': rust={rust_d} cpp={cpp_d}" + ); + } +} + +#[test] +fn ffi_spearman_dist_profile_equivalence() { + let ctx = DistanceContext::new(4); + let cpp = OriginalSpearmanDist::new(4); + + let rust_kmap = KmerMap::new(4, true); + let seqs = ["ACGTACGTACGTACGT", "AAAAAAAAAAAAAAAA", "ATCGATCGATCGATCG"]; + + let profiles: Vec> = seqs + .iter() + .map(|s| { + Profiler::new(4, s.as_bytes(), &rust_kmap) + .get_profile() + .to_vec() + }) + .collect(); + + for i in 0..profiles.len() { + for j in (i + 1)..profiles.len() { + let rust_d = spearman_dist_profiles(&ctx, &profiles[i], &profiles[j]); + let cpp_d = cpp.get_dist_profile(&profiles[i], &profiles[j]); + + let diff = (rust_d - cpp_d).abs(); + assert!( + diff < 1e-10, + "SpearmanDist profile mismatch: rust={rust_d} cpp={cpp_d}" + ); + } + } +} + +#[test] +fn ffi_spearman_normalized_equivalence() { + let mut ctx = DistanceContext::new(4); + ctx.set_normalization(true); + let cpp = OriginalSpearmanDist::new(4); + cpp.set_normalization(true); + + let seq1 = "ACGTACGTACGTACGT"; + let seq2 = "AAAAAAAAAAAAAAAA"; + + let rust_d = spearman_dist_seq(&ctx, seq1.as_bytes(), seq2.as_bytes()); + let cpp_d = cpp.get_dist_seq(seq1, seq2); + + let diff = (rust_d - cpp_d).abs(); + assert!( + diff < 1e-10, + "SpearmanDist normalized mismatch: rust={rust_d} cpp={cpp_d}" + ); +} diff --git a/crates/maxbin-rs/tests/emanager_equivalence.rs b/crates/maxbin-rs-equivalence/tests/emanager_equivalence.rs similarity index 96% rename from crates/maxbin-rs/tests/emanager_equivalence.rs rename to crates/maxbin-rs-equivalence/tests/emanager_equivalence.rs index 806d160..a766b29 100644 --- a/crates/maxbin-rs/tests/emanager_equivalence.rs +++ b/crates/maxbin-rs-equivalence/tests/emanager_equivalence.rs @@ -61,8 +61,11 @@ fn emanager_cpp_vs_rust() { // Run C++ EManager via FFI. let cpp_prefix = test_dir.join("cpp").to_str().unwrap().to_string(); { - let em = - maxbin_rs::original_ffi::OriginalEManager::new(&fasta_path, &abund_path, &cpp_prefix); + let em = maxbin_rs_equivalence::original_ffi::OriginalEManager::new( + &fasta_path, + &abund_path, + &cpp_prefix, + ); em.set_thread_num(1); assert_eq!(em.run(&seed_path), 0, "C++ EManager returned error"); } @@ -140,8 +143,11 @@ fn emanager_precision_divergence_cami_depth4() { // Run C++ EManager via FFI let cpp_prefix = test_dir.join("cpp").to_str().unwrap().to_string(); { - let em = - maxbin_rs::original_ffi::OriginalEManager::new(&fasta_path, &abund_path, &cpp_prefix); + let em = maxbin_rs_equivalence::original_ffi::OriginalEManager::new( + &fasta_path, + &abund_path, + &cpp_prefix, + ); em.set_thread_num(1); let rc = em.run(&seed_path); assert_eq!(rc, 0, "C++ EManager returned error"); diff --git a/crates/maxbin-rs-equivalence/tests/fasta_ffi.rs b/crates/maxbin-rs-equivalence/tests/fasta_ffi.rs new file mode 100644 index 0000000..e36041c --- /dev/null +++ b/crates/maxbin-rs-equivalence/tests/fasta_ffi.rs @@ -0,0 +1,71 @@ +//! FFI equivalence test for the FASTA parser. +//! +//! Compares `maxbin_rs::fasta` against the original MaxBin2 C++ fastaReader. +//! Requires nf-core test data via the `MAXBIN2_TEST_CONTIGS` env var (set by +//! the Nix devshell); skips otherwise. Previously an in-source test in +//! `fasta.rs`. +use maxbin_rs::fasta::parse_file; +use maxbin_rs_equivalence::original_ffi::OriginalFastaReader; + +/// Decompress a .gz file to a temp path for the C++ parser (which can't read gzip). +fn decompress_gz(gz_path: &std::path::Path) -> std::path::PathBuf { + use std::io::{Read, Write}; + let out_path = std::env::temp_dir().join("maxbin_rs_test_contigs.fa"); + let file = std::fs::File::open(gz_path).unwrap(); + let mut decoder = flate2::read::GzDecoder::new(file); + let mut contents = Vec::new(); + decoder.read_to_end(&mut contents).unwrap(); + let mut out = std::fs::File::create(&out_path).unwrap(); + out.write_all(&contents).unwrap(); + out_path +} + +#[test] +fn equivalence_with_original_fasta_parser() { + // This test requires nf-core test data (only available in devshell). + let gz_path = match std::env::var("MAXBIN2_TEST_CONTIGS") { + Ok(p) => std::path::PathBuf::from(p), + Err(_) => return, + }; + + // Parse with our Rust implementation (handles gzip natively) + let rust_records = parse_file(&gz_path).expect("rust parser failed"); + + // Decompress for the C++ parser, then parse with the original + let plain_path = decompress_gz(&gz_path); + let original = OriginalFastaReader::new(&plain_path); + + // Same number of records + assert_eq!( + rust_records.len(), + original.num_records() as usize, + "record count mismatch" + ); + + // Field-by-field comparison for every record + for (i, rust_rec) in rust_records.iter().enumerate() { + let idx = i as u32; + + assert_eq!( + rust_rec.header, + original.header(idx), + "header mismatch at record {i}" + ); + + assert_eq!( + rust_rec.seq, + original.seq(idx), + "sequence mismatch at record {i} (header: {})", + rust_rec.header + ); + + assert_eq!( + rust_rec.len(), + original.seq_len(idx) as usize, + "length mismatch at record {i} (header: {})", + rust_rec.header + ); + } + + std::fs::remove_file(&plain_path).ok(); +} diff --git a/crates/maxbin-rs/tests/fixtures/divergent-em/abund b/crates/maxbin-rs-equivalence/tests/fixtures/divergent-em/abund similarity index 100% rename from crates/maxbin-rs/tests/fixtures/divergent-em/abund rename to crates/maxbin-rs-equivalence/tests/fixtures/divergent-em/abund diff --git a/crates/maxbin-rs/tests/fixtures/divergent-em/contigs.fa b/crates/maxbin-rs-equivalence/tests/fixtures/divergent-em/contigs.fa similarity index 100% rename from crates/maxbin-rs/tests/fixtures/divergent-em/contigs.fa rename to crates/maxbin-rs-equivalence/tests/fixtures/divergent-em/contigs.fa diff --git a/crates/maxbin-rs/tests/fixtures/divergent-em/seed b/crates/maxbin-rs-equivalence/tests/fixtures/divergent-em/seed similarity index 100% rename from crates/maxbin-rs/tests/fixtures/divergent-em/seed rename to crates/maxbin-rs-equivalence/tests/fixtures/divergent-em/seed diff --git a/crates/maxbin-rs-equivalence/tests/kmer_map_ffi.rs b/crates/maxbin-rs-equivalence/tests/kmer_map_ffi.rs new file mode 100644 index 0000000..5c4dc08 --- /dev/null +++ b/crates/maxbin-rs-equivalence/tests/kmer_map_ffi.rs @@ -0,0 +1,55 @@ +//! FFI equivalence tests for the k-mer mapping table. +//! +//! Compares `maxbin_rs::kmer_map::KmerMap` against the original MaxBin2 C++ +//! kmerMap. Previously in-source tests in `kmer_map.rs`. +use maxbin_rs::kmer_map::KmerMap; +use maxbin_rs_equivalence::original_ffi::OriginalKmerMap; + +#[test] +fn ffi_equivalence_kmer4_symmetric() { + let rust_km = KmerMap::new(4, true); + let cpp_km = OriginalKmerMap::new(4, true); + + assert_eq!(rust_km.get_entry_num() as i32, cpp_km.entry_num()); + + let cpp_table = cpp_km.get_full_table(); + let rust_table = rust_km.get_full_table(); + + assert_eq!(rust_table.len(), cpp_table.len()); + for i in 0..rust_table.len() { + assert_eq!(rust_table[i], cpp_table[i], "mapping mismatch at index {i}"); + } +} + +#[test] +fn ffi_equivalence_kmer4_nonsymmetric() { + let rust_km = KmerMap::new(4, false); + let cpp_km = OriginalKmerMap::new(4, false); + + assert_eq!(rust_km.get_entry_num() as i32, cpp_km.entry_num()); + + let cpp_table = cpp_km.get_full_table(); + let rust_table = rust_km.get_full_table(); + + for i in 0..rust_table.len() { + assert_eq!(rust_table[i], cpp_table[i]); + } +} + +#[test] +fn ffi_equivalence_kmer3_symmetric() { + let rust_km = KmerMap::new(3, true); + let cpp_km = OriginalKmerMap::new(3, true); + + assert_eq!(rust_km.get_entry_num() as i32, cpp_km.entry_num()); + + let cpp_table = cpp_km.get_full_table(); + let rust_table = rust_km.get_full_table(); + + for i in 0..rust_table.len() { + assert_eq!( + rust_table[i], cpp_table[i], + "kmer3 symmetric mapping mismatch at index {i}" + ); + } +} diff --git a/crates/maxbin-rs-equivalence/tests/normal_distribution_ffi.rs b/crates/maxbin-rs-equivalence/tests/normal_distribution_ffi.rs new file mode 100644 index 0000000..d53c3a3 --- /dev/null +++ b/crates/maxbin-rs-equivalence/tests/normal_distribution_ffi.rs @@ -0,0 +1,40 @@ +//! FFI equivalence test for the normal distribution PDF. +//! +//! Compares `maxbin_rs::normal_distribution::NormalDistribution` against the +//! original MaxBin2 C++ NormalDistribution. Previously an in-source test in +//! `normal_distribution.rs`. +use maxbin_rs::normal_distribution::NormalDistribution; +use maxbin_rs_equivalence::original_ffi::OriginalNormalDistribution; + +#[test] +fn ffi_equivalence() { + let test_cases = [ + (0.0, 1.0, 0.0), + (0.0, 1.0, 1.0), + (0.0, 1.0, -1.0), + (5.0, 2.0, 5.0), + (5.0, 2.0, 3.0), + (5.0, 2.0, 7.0), + (10.0, 0.5, 10.0), + (10.0, 0.5, 9.5), + (-3.0, 1.5, -3.0), + (-3.0, 1.5, 0.0), + (0.0, 1.0, 3.5), + (100.0, 10.0, 95.0), + ]; + + for &(mean, std, input) in &test_cases { + let rust_nd = NormalDistribution::new(mean, std); + let cpp_nd = OriginalNormalDistribution::new(mean, std); + + let rust_p = rust_nd.prob(input); + let cpp_p = cpp_nd.prob(input); + + assert_eq!( + rust_p.to_bits(), + cpp_p.to_bits(), + "NormalDist not bit-identical for mean={mean}, std={std}, input={input}: \ + rust={rust_p:e} cpp={cpp_p:e}" + ); + } +} diff --git a/crates/maxbin-rs-equivalence/tests/profiler_ffi.rs b/crates/maxbin-rs-equivalence/tests/profiler_ffi.rs new file mode 100644 index 0000000..0bc2fd3 --- /dev/null +++ b/crates/maxbin-rs-equivalence/tests/profiler_ffi.rs @@ -0,0 +1,48 @@ +//! FFI equivalence test for the tetranucleotide profiler. +//! +//! Compares `maxbin_rs::profiler::Profiler` against the original MaxBin2 C++ +//! Profiler. Previously an in-source test in `profiler.rs`. +use maxbin_rs::kmer_map::KmerMap; +use maxbin_rs::profiler::Profiler; +use maxbin_rs_equivalence::original_ffi::{OriginalKmerMap, OriginalProfiler}; + +#[test] +fn ffi_equivalence() { + let test_seqs = [ + "ACGTACGTACGTACGT", + "AAAAAAAAAAAA", + "ACGTTTTTGGGGCCCC", + "ATCGATCGATCGATCGATCG", + "ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGT", + ]; + + let rust_kmap = KmerMap::new(4, true); + let cpp_kmap = OriginalKmerMap::new(4, true); + let entry_num = rust_kmap.get_entry_num(); + + for seq in &test_seqs { + let rust_prof = Profiler::new(4, seq.as_bytes(), &rust_kmap); + let cpp_prof = OriginalProfiler::new(4, seq, &cpp_kmap); + + let rust_profile = rust_prof.get_profile(); + let cpp_profile = cpp_prof.get_profile(entry_num as i32); + + assert_eq!(rust_profile.len(), cpp_profile.len()); + for i in 0..rust_profile.len() { + assert_eq!( + rust_profile[i].to_bits(), + cpp_profile[i].to_bits(), + "profile not bit-identical at index {i} for seq '{seq}': rust={:e} cpp={:e}", + rust_profile[i], + cpp_profile[i] + ); + } + + let rust_pn = rust_prof.get_percent_n(); + let cpp_pn = cpp_prof.get_percent_n(); + assert!( + (rust_pn - cpp_pn).abs() < 1e-6, + "percent_n mismatch for seq '{seq}': rust={rust_pn} cpp={cpp_pn}" + ); + } +} diff --git a/crates/maxbin-rs/tests/proptest_abundance.rs b/crates/maxbin-rs-equivalence/tests/proptest_abundance.rs similarity index 97% rename from crates/maxbin-rs/tests/proptest_abundance.rs rename to crates/maxbin-rs-equivalence/tests/proptest_abundance.rs index 016c1ac..29794b9 100644 --- a/crates/maxbin-rs/tests/proptest_abundance.rs +++ b/crates/maxbin-rs-equivalence/tests/proptest_abundance.rs @@ -70,7 +70,7 @@ proptest! { let rust_records = maxbin_rs::abundance::parse(Cursor::new(input.as_bytes())).unwrap(); // C++ - let original = maxbin_rs::original_ffi::OriginalAbundanceLoader::new(&tmp_path); + let original = maxbin_rs_equivalence::original_ffi::OriginalAbundanceLoader::new(&tmp_path); prop_assert!(original.is_parse_success(), "C++ parse failed"); prop_assert_eq!( diff --git a/crates/maxbin-rs/tests/proptest_distance.rs b/crates/maxbin-rs-equivalence/tests/proptest_distance.rs similarity index 95% rename from crates/maxbin-rs/tests/proptest_distance.rs rename to crates/maxbin-rs-equivalence/tests/proptest_distance.rs index df0e51a..7fd1d19 100644 --- a/crates/maxbin-rs/tests/proptest_distance.rs +++ b/crates/maxbin-rs-equivalence/tests/proptest_distance.rs @@ -41,7 +41,7 @@ mod distance_proptest { #[test] fn euc_dist_equivalence(seq1 in arb_dna_seq(), seq2 in arb_dna_seq()) { let ctx = maxbin_rs::distance::DistanceContext::new(4); - let cpp = maxbin_rs::original_ffi::OriginalEucDist::new(4); + let cpp = maxbin_rs_equivalence::original_ffi::OriginalEucDist::new(4); let rust_d = maxbin_rs::distance::euc_dist_seq(&ctx, seq1.as_bytes(), seq2.as_bytes()); let cpp_d = cpp.get_dist_seq(&seq1, &seq2); @@ -57,7 +57,7 @@ mod distance_proptest { #[test] fn spearman_dist_equivalence(seq1 in arb_dna_seq(), seq2 in arb_dna_seq()) { let ctx = maxbin_rs::distance::DistanceContext::new(4); - let cpp = maxbin_rs::original_ffi::OriginalSpearmanDist::new(4); + let cpp = maxbin_rs_equivalence::original_ffi::OriginalSpearmanDist::new(4); let rust_d = maxbin_rs::distance::spearman_dist_seq(&ctx, seq1.as_bytes(), seq2.as_bytes()); let cpp_d = cpp.get_dist_seq(&seq1, &seq2); @@ -125,7 +125,7 @@ mod distance_profiles_proptest { let p2 = prof2.get_profile().to_vec(); let rust_d = maxbin_rs::distance::euc_dist_profiles(&p1, &p2); - let cpp = maxbin_rs::original_ffi::OriginalEucDist::new(4); + let cpp = maxbin_rs_equivalence::original_ffi::OriginalEucDist::new(4); let cpp_d = cpp.get_dist_profile(&p1, &p2); let diff = (rust_d - cpp_d).abs(); @@ -148,7 +148,7 @@ mod distance_profiles_proptest { let ctx = maxbin_rs::distance::DistanceContext::new(4); let rust_d = maxbin_rs::distance::spearman_dist_profiles(&ctx, &p1, &p2); - let cpp = maxbin_rs::original_ffi::OriginalSpearmanDist::new(4); + let cpp = maxbin_rs_equivalence::original_ffi::OriginalSpearmanDist::new(4); let cpp_d = cpp.get_dist_profile(&p1, &p2); let diff = (rust_d - cpp_d).abs(); diff --git a/crates/maxbin-rs/tests/proptest_emanager.rs b/crates/maxbin-rs-equivalence/tests/proptest_emanager.rs similarity index 95% rename from crates/maxbin-rs/tests/proptest_emanager.rs rename to crates/maxbin-rs-equivalence/tests/proptest_emanager.rs index 5766a75..7116b8e 100644 --- a/crates/maxbin-rs/tests/proptest_emanager.rs +++ b/crates/maxbin-rs-equivalence/tests/proptest_emanager.rs @@ -15,8 +15,10 @@ use proptest::prelude::*; /// Compute prob_dist via C++ FFI NormalDistribution for cross-checking. fn prob_dist_via_cpp_ffi(distance: f64) -> f64 { - let intra = maxbin_rs::original_ffi::OriginalNormalDistribution::new(0.0, 0.01037897 / 2.0); - let inter = maxbin_rs::original_ffi::OriginalNormalDistribution::new(0.0676654, 0.03419337); + let intra = + maxbin_rs_equivalence::original_ffi::OriginalNormalDistribution::new(0.0, 0.01037897 / 2.0); + let inter = + maxbin_rs_equivalence::original_ffi::OriginalNormalDistribution::new(0.0676654, 0.03419337); let d_intra = intra.prob(distance); let d_inter = inter.prob(distance); d_intra / (d_inter + d_intra) diff --git a/crates/maxbin-rs/tests/proptest_fasta.rs b/crates/maxbin-rs-equivalence/tests/proptest_fasta.rs similarity index 97% rename from crates/maxbin-rs/tests/proptest_fasta.rs rename to crates/maxbin-rs-equivalence/tests/proptest_fasta.rs index 4dab953..7320a93 100644 --- a/crates/maxbin-rs/tests/proptest_fasta.rs +++ b/crates/maxbin-rs-equivalence/tests/proptest_fasta.rs @@ -79,7 +79,7 @@ proptest! { let rust_records = maxbin_rs::fasta::parse(Cursor::new(input.as_bytes())); // C++ parser - let original = maxbin_rs::original_ffi::OriginalFastaReader::new(&tmp_path); + let original = maxbin_rs_equivalence::original_ffi::OriginalFastaReader::new(&tmp_path); prop_assert_eq!( rust_records.len(), diff --git a/crates/maxbin-rs/tests/proptest_kmermap.rs b/crates/maxbin-rs-equivalence/tests/proptest_kmermap.rs similarity index 91% rename from crates/maxbin-rs/tests/proptest_kmermap.rs rename to crates/maxbin-rs-equivalence/tests/proptest_kmermap.rs index e5bfc58..47fadc2 100644 --- a/crates/maxbin-rs/tests/proptest_kmermap.rs +++ b/crates/maxbin-rs-equivalence/tests/proptest_kmermap.rs @@ -14,7 +14,7 @@ mod kmermap_proptest { #[test] fn kmermap_mapping_equivalence(kmerlen in 1..6i32, symmetric in proptest::bool::ANY) { let rust_km = maxbin_rs::kmer_map::KmerMap::new(kmerlen as usize, symmetric); - let cpp_km = maxbin_rs::original_ffi::OriginalKmerMap::new(kmerlen, symmetric); + let cpp_km = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(kmerlen, symmetric); let rust_entry = rust_km.get_entry_num() as i32; let cpp_entry = cpp_km.entry_num(); @@ -64,7 +64,7 @@ mod kmermap_per_kmer_proptest { fn kmermap_get_mapping(args in arb_kmerlen_and_kmer()) { let (kmerlen, symmetric, kmer) = args; let rust_km = maxbin_rs::kmer_map::KmerMap::new(kmerlen, symmetric); - let cpp_km = maxbin_rs::original_ffi::OriginalKmerMap::new(kmerlen as i32, symmetric); + let cpp_km = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(kmerlen as i32, symmetric); let rust_idx = rust_km.get_mapping(kmer.as_bytes()); let cpp_idx = cpp_km.get_mapping(&kmer); @@ -81,7 +81,7 @@ mod kmermap_per_kmer_proptest { fn kmermap_get_reverse_mapping_str(args in arb_kmerlen_and_kmer()) { let (kmerlen, symmetric, kmer) = args; let rust_km = maxbin_rs::kmer_map::KmerMap::new(kmerlen, symmetric); - let cpp_km = maxbin_rs::original_ffi::OriginalKmerMap::new(kmerlen as i32, symmetric); + let cpp_km = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(kmerlen as i32, symmetric); let rust_idx = rust_km.get_reverse_mapping_str(kmer.as_bytes()); let cpp_idx = cpp_km.get_reverse_mapping_str(&kmer); @@ -97,7 +97,7 @@ mod kmermap_per_kmer_proptest { fn kmermap_get_reverse_mapping_idx(args in arb_kmerlen_and_index()) { let (kmerlen, symmetric, index) = args; let rust_km = maxbin_rs::kmer_map::KmerMap::new(kmerlen, symmetric); - let cpp_km = maxbin_rs::original_ffi::OriginalKmerMap::new(kmerlen as i32, symmetric); + let cpp_km = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(kmerlen as i32, symmetric); let rust_result = rust_km.get_reverse_mapping_idx(index); let cpp_result = cpp_km.get_reverse_mapping_idx(index as i32); @@ -115,7 +115,7 @@ mod kmermap_per_kmer_proptest { kmer in "[A-Za-z0-9NRYWSKM]{4}" ) { let rust_km = maxbin_rs::kmer_map::KmerMap::new(4, symmetric); - let cpp_km = maxbin_rs::original_ffi::OriginalKmerMap::new(4, symmetric); + let cpp_km = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(4, symmetric); let rust_idx = rust_km.get_mapping(kmer.as_bytes()); let cpp_idx = cpp_km.get_mapping(&kmer); @@ -133,7 +133,7 @@ mod kmermap_per_kmer_proptest { kmer in "[acgtACGT]{4}" ) { let rust_km = maxbin_rs::kmer_map::KmerMap::new(4, symmetric); - let cpp_km = maxbin_rs::original_ffi::OriginalKmerMap::new(4, symmetric); + let cpp_km = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(4, symmetric); let rust_idx = rust_km.get_mapping(kmer.as_bytes()); let cpp_idx = cpp_km.get_mapping(&kmer); diff --git a/crates/maxbin-rs/tests/proptest_normal_distribution.rs b/crates/maxbin-rs-equivalence/tests/proptest_normal_distribution.rs similarity index 91% rename from crates/maxbin-rs/tests/proptest_normal_distribution.rs rename to crates/maxbin-rs-equivalence/tests/proptest_normal_distribution.rs index 073cf19..bcc9000 100644 --- a/crates/maxbin-rs/tests/proptest_normal_distribution.rs +++ b/crates/maxbin-rs-equivalence/tests/proptest_normal_distribution.rs @@ -14,7 +14,7 @@ proptest! { input in -200.0f64..200.0 ) { let rust_nd = maxbin_rs::normal_distribution::NormalDistribution::new(mean, std_dev); - let cpp_nd = maxbin_rs::original_ffi::OriginalNormalDistribution::new(mean, std_dev); + let cpp_nd = maxbin_rs_equivalence::original_ffi::OriginalNormalDistribution::new(mean, std_dev); let rust_p = rust_nd.prob(input); let cpp_p = cpp_nd.prob(input); diff --git a/crates/maxbin-rs/tests/proptest_profiler.rs b/crates/maxbin-rs-equivalence/tests/proptest_profiler.rs similarity index 98% rename from crates/maxbin-rs/tests/proptest_profiler.rs rename to crates/maxbin-rs-equivalence/tests/proptest_profiler.rs index 4793512..47c2ad5 100644 --- a/crates/maxbin-rs/tests/proptest_profiler.rs +++ b/crates/maxbin-rs-equivalence/tests/proptest_profiler.rs @@ -48,11 +48,11 @@ mod profiler_proptest { #[test] fn profiler_equivalence(seq in arb_dna_seq()) { let rust_kmap = maxbin_rs::kmer_map::KmerMap::new(4, true); - let cpp_kmap = maxbin_rs::original_ffi::OriginalKmerMap::new(4, true); + let cpp_kmap = maxbin_rs_equivalence::original_ffi::OriginalKmerMap::new(4, true); let entry_num = rust_kmap.get_entry_num(); let rust_prof = maxbin_rs::profiler::Profiler::new(4, seq.as_bytes(), &rust_kmap); - let cpp_prof = maxbin_rs::original_ffi::OriginalProfiler::new(4, &seq, &cpp_kmap); + let cpp_prof = maxbin_rs_equivalence::original_ffi::OriginalProfiler::new(4, &seq, &cpp_kmap); let rust_profile = rust_prof.get_profile(); let cpp_profile = cpp_prof.get_profile(entry_num as i32); diff --git a/crates/maxbin-rs/tests/proptest_quicksort.rs b/crates/maxbin-rs-equivalence/tests/proptest_quicksort.rs similarity index 88% rename from crates/maxbin-rs/tests/proptest_quicksort.rs rename to crates/maxbin-rs-equivalence/tests/proptest_quicksort.rs index 1565a00..ee866b9 100644 --- a/crates/maxbin-rs/tests/proptest_quicksort.rs +++ b/crates/maxbin-rs-equivalence/tests/proptest_quicksort.rs @@ -17,7 +17,7 @@ proptest! { // C++ let mut orig_arr = input.clone(); let mut orig_idx: Vec = (0..input.len() as i32).collect(); - maxbin_rs::original_ffi::original_quicksort(&mut orig_arr, Some(&mut orig_idx)); + maxbin_rs_equivalence::original_ffi::original_quicksort(&mut orig_arr, Some(&mut orig_idx)); prop_assert_eq!(&rust_arr, &orig_arr, "array mismatch"); prop_assert_eq!(&rust_idx, &orig_idx, "index mismatch"); @@ -30,7 +30,7 @@ proptest! { maxbin_rs::quicksort::sort_descending(&mut rust_arr, None); let mut orig_arr = input.clone(); - maxbin_rs::original_ffi::original_quicksort(&mut orig_arr, None); + maxbin_rs_equivalence::original_ffi::original_quicksort(&mut orig_arr, None); prop_assert_eq!(&rust_arr, &orig_arr, "array mismatch (no indices)"); } @@ -49,7 +49,7 @@ proptest! { let mut orig_arr = input.clone(); let mut orig_idx: Vec = (0..input.len() as i32).collect(); - maxbin_rs::original_ffi::original_quicksort(&mut orig_arr, Some(&mut orig_idx)); + maxbin_rs_equivalence::original_ffi::original_quicksort(&mut orig_arr, Some(&mut orig_idx)); prop_assert_eq!(&rust_arr, &orig_arr, "array mismatch with duplicates"); prop_assert_eq!(&rust_idx, &orig_idx, "index mismatch with duplicates"); diff --git a/crates/maxbin-rs-equivalence/tests/quicksort_ffi.rs b/crates/maxbin-rs-equivalence/tests/quicksort_ffi.rs new file mode 100644 index 0000000..e18571f --- /dev/null +++ b/crates/maxbin-rs-equivalence/tests/quicksort_ffi.rs @@ -0,0 +1,43 @@ +//! FFI equivalence test for the descending quicksort. +//! +//! Compares `maxbin_rs::quicksort::sort_descending` against the original +//! MaxBin2 C++ quickSort. Previously an in-source test in `quicksort.rs`. +use maxbin_rs::quicksort::sort_descending; +use maxbin_rs_equivalence::original_ffi::original_quicksort; + +#[test] +fn equivalence_with_original() { + // Test a variety of inputs against the C++ quickSort via FFI. + let test_cases: Vec> = vec![ + vec![1.0, 5.0, 3.0, 2.0, 4.0], + vec![5.0, 4.0, 3.0, 2.0, 1.0], + vec![1.0, 1.0, 1.0], + vec![3.0, 1.0, 4.0, 1.0, 5.0, 9.0, 2.0, 6.0, 5.0, 3.0, 5.0], + vec![42.0], + vec![1.0, 2.0], + vec![2.0, 1.0], + vec![0.001, 0.01, 0.1, 1.0, 10.0, 100.0], + vec![-1.0, 0.0, 1.0, -0.5, 0.5], + ]; + + for (case_idx, case) in test_cases.iter().enumerate() { + // Rust version + let mut rust_arr = case.clone(); + let mut rust_idx: Vec = (0..case.len() as i32).collect(); + sort_descending(&mut rust_arr, Some(&mut rust_idx)); + + // Original C++ version + let mut orig_arr = case.clone(); + let mut orig_idx: Vec = (0..case.len() as i32).collect(); + original_quicksort(&mut orig_arr, Some(&mut orig_idx)); + + assert_eq!( + rust_arr, orig_arr, + "array mismatch for case {case_idx}: {case:?}" + ); + assert_eq!( + rust_idx, orig_idx, + "index mismatch for case {case_idx}: {case:?}" + ); + } +} diff --git a/crates/maxbin-rs/vendor/ffi/abundance_ffi.cpp b/crates/maxbin-rs-equivalence/vendor/ffi/abundance_ffi.cpp similarity index 100% rename from crates/maxbin-rs/vendor/ffi/abundance_ffi.cpp rename to crates/maxbin-rs-equivalence/vendor/ffi/abundance_ffi.cpp diff --git a/crates/maxbin-rs/vendor/ffi/distance_ffi.cpp b/crates/maxbin-rs-equivalence/vendor/ffi/distance_ffi.cpp similarity index 100% rename from crates/maxbin-rs/vendor/ffi/distance_ffi.cpp rename to crates/maxbin-rs-equivalence/vendor/ffi/distance_ffi.cpp diff --git a/crates/maxbin-rs/vendor/ffi/emanager_ffi.cpp b/crates/maxbin-rs-equivalence/vendor/ffi/emanager_ffi.cpp similarity index 100% rename from crates/maxbin-rs/vendor/ffi/emanager_ffi.cpp rename to crates/maxbin-rs-equivalence/vendor/ffi/emanager_ffi.cpp diff --git a/crates/maxbin-rs/vendor/ffi/fasta_ffi.cpp b/crates/maxbin-rs-equivalence/vendor/ffi/fasta_ffi.cpp similarity index 100% rename from crates/maxbin-rs/vendor/ffi/fasta_ffi.cpp rename to crates/maxbin-rs-equivalence/vendor/ffi/fasta_ffi.cpp diff --git a/crates/maxbin-rs/vendor/ffi/kmermap_ffi.cpp b/crates/maxbin-rs-equivalence/vendor/ffi/kmermap_ffi.cpp similarity index 100% rename from crates/maxbin-rs/vendor/ffi/kmermap_ffi.cpp rename to crates/maxbin-rs-equivalence/vendor/ffi/kmermap_ffi.cpp diff --git a/crates/maxbin-rs/vendor/ffi/normaldist_ffi.cpp b/crates/maxbin-rs-equivalence/vendor/ffi/normaldist_ffi.cpp similarity index 100% rename from crates/maxbin-rs/vendor/ffi/normaldist_ffi.cpp rename to crates/maxbin-rs-equivalence/vendor/ffi/normaldist_ffi.cpp diff --git a/crates/maxbin-rs/vendor/ffi/profiler_ffi.cpp b/crates/maxbin-rs-equivalence/vendor/ffi/profiler_ffi.cpp similarity index 100% rename from crates/maxbin-rs/vendor/ffi/profiler_ffi.cpp rename to crates/maxbin-rs-equivalence/vendor/ffi/profiler_ffi.cpp diff --git a/crates/maxbin-rs/vendor/ffi/quicksort_ffi.cpp b/crates/maxbin-rs-equivalence/vendor/ffi/quicksort_ffi.cpp similarity index 100% rename from crates/maxbin-rs/vendor/ffi/quicksort_ffi.cpp rename to crates/maxbin-rs-equivalence/vendor/ffi/quicksort_ffi.cpp diff --git a/crates/maxbin-rs/vendor/maxbin2-cpp-ffi.patch b/crates/maxbin-rs-equivalence/vendor/maxbin2-cpp-ffi.patch similarity index 100% rename from crates/maxbin-rs/vendor/maxbin2-cpp-ffi.patch rename to crates/maxbin-rs-equivalence/vendor/maxbin2-cpp-ffi.patch diff --git a/crates/maxbin-rs/Cargo.toml b/crates/maxbin-rs/Cargo.toml index 8e8fb32..05e9c2d 100644 --- a/crates/maxbin-rs/Cargo.toml +++ b/crates/maxbin-rs/Cargo.toml @@ -12,9 +12,3 @@ serde = { version = "1", features = ["derive"] } serde_json = "1" seq_io = "0.3" tempfile = "3" - -[build-dependencies] -cc = "1.2.60" - -[dev-dependencies] -proptest = "1.11.0" diff --git a/crates/maxbin-rs/src/abundance.rs b/crates/maxbin-rs/src/abundance.rs index 96817cd..d897d8f 100644 --- a/crates/maxbin-rs/src/abundance.rs +++ b/crates/maxbin-rs/src/abundance.rs @@ -152,41 +152,4 @@ mod tests { assert_eq!(records[1].header, "c2"); assert_eq!(records[1].abundance, 2.5); } - - #[test] - fn equivalence_with_original() { - // Writes a temp abundance file, parses with both Rust and C++ FFI, - // compares field by field. - use std::io::Write; - - let test_data = "contig_1\t10.5\ncontig_2\t0.00001\ncontig_3\t100.0\ncontig_4 0.5\n"; - - let tmp_path = std::env::temp_dir().join("maxbin_rs_test_abund.txt"); - { - let mut f = std::fs::File::create(&tmp_path).unwrap(); - f.write_all(test_data.as_bytes()).unwrap(); - } - - // Rust - let rust_records = parse_file(&tmp_path).unwrap(); - - // Original C++ - let original = crate::original_ffi::OriginalAbundanceLoader::new(&tmp_path); - - assert_eq!(rust_records.len(), original.num_records() as usize); - assert!(original.is_parse_success()); - - for (i, rec) in rust_records.iter().enumerate() { - let orig_abund = original.abundance_by_index(i as i32); - assert!( - (rec.abundance - orig_abund).abs() < 1e-10, - "abundance mismatch at record {i} (header: {}): rust={} original={}", - rec.header, - rec.abundance, - orig_abund, - ); - } - - std::fs::remove_file(&tmp_path).ok(); - } } diff --git a/crates/maxbin-rs/src/cli.rs b/crates/maxbin-rs/src/cli.rs index e19f264..79c7938 100644 --- a/crates/maxbin-rs/src/cli.rs +++ b/crates/maxbin-rs/src/cli.rs @@ -9,7 +9,6 @@ /// maxbin-rs filter — filter contigs by minimum length /// maxbin-rs seeds — generate seed file from HMMER marker gene hits /// maxbin-rs em — run the Rust EM algorithm -/// maxbin-rs cpp-em — run the original C++ EM via FFI (equivalence testing) /// maxbin-rs sam-to-abund — compute abundance from a SAM file /// maxbin-rs pipeline — run the full pipeline (default when no subcommand given) use std::path::PathBuf; @@ -39,9 +38,6 @@ pub enum Command { Seeds(SeedsArgs), /// Run the Rust EM algorithm on pre-computed inputs. Em(EmArgs), - /// Run the original C++ EM via FFI (for equivalence testing). - #[command(name = "cpp-em", alias = "cpp_em")] - CppEm(CppEmArgs), /// Compute abundance from a SAM file. #[command(name = "sam-to-abund", alias = "sam_to_abund")] SamToAbund(SamToAbundArgs), @@ -159,20 +155,6 @@ pub struct PipelineArgs { pub work_dir: Option, } -#[derive(Debug, Args)] -pub struct CppEmArgs { - #[arg(long)] - pub contig: PathBuf, - #[arg(long)] - pub abund: PathBuf, - #[arg(long)] - pub seed: PathBuf, - #[arg(long)] - pub out: String, - #[arg(long, default_value_t = 1)] - pub thread: usize, -} - #[derive(Debug, Args)] pub struct SamToAbundArgs { #[arg(long)] diff --git a/crates/maxbin-rs/src/distance.rs b/crates/maxbin-rs/src/distance.rs index 5a29ab4..589166e 100644 --- a/crates/maxbin-rs/src/distance.rs +++ b/crates/maxbin-rs/src/distance.rs @@ -192,138 +192,4 @@ mod tests { assert!(d >= 0.0); assert!(d <= 1.0, "normalized distance should be <= 1, got {d}"); } - - // --- FFI equivalence tests --- - - #[test] - fn ffi_euc_dist_seq_equivalence() { - let test_pairs = [ - ("ACGTACGTACGTACGT", "ACGTACGTACGTACGT"), - ("AAAAAAAAAAAAAAAA", "CCCCCCCCCCCCCCCC"), - ("ATCGATCGATCGATCG", "TAGCTAGCTAGCTAGC"), - ("ACGTACGTACGTACGTACGTACGT", "TTTTAAAACCCCGGGG"), - ("AACCGGTTAACCGGTT", "GGTTAACCGGTTAACC"), - ]; - - let ctx = DistanceContext::new(4); - let cpp = crate::original_ffi::OriginalEucDist::new(4); - - for (seq1, seq2) in &test_pairs { - let rust_d = euc_dist_seq(&ctx, seq1.as_bytes(), seq2.as_bytes()); - let cpp_d = cpp.get_dist_seq(seq1, seq2); - - assert_eq!( - rust_d.to_bits(), - cpp_d.to_bits(), - "EucDist seq not bit-identical for '{seq1}' vs '{seq2}': rust={rust_d:e} cpp={cpp_d:e}" - ); - } - } - - #[test] - fn ffi_euc_dist_profile_equivalence() { - // Build profiles from known sequences, then compare profile-based distance - let cpp = crate::original_ffi::OriginalEucDist::new(4); - - let seqs = ["ACGTACGTACGTACGT", "AAAAAAAAAAAAAAAA", "ATCGATCGATCGATCG"]; - - let rust_kmap = KmerMap::new(4, true); - let profiles: Vec> = seqs - .iter() - .map(|s| { - Profiler::new(4, s.as_bytes(), &rust_kmap) - .get_profile() - .to_vec() - }) - .collect(); - - for i in 0..profiles.len() { - for j in (i + 1)..profiles.len() { - let rust_d = euc_dist_profiles(&profiles[i], &profiles[j]); - let cpp_d = cpp.get_dist_profile(&profiles[i], &profiles[j]); - - assert_eq!( - rust_d.to_bits(), - cpp_d.to_bits(), - "EucDist profile not bit-identical: rust={rust_d:e} cpp={cpp_d:e}" - ); - } - } - } - - #[test] - fn ffi_spearman_dist_seq_equivalence() { - let test_pairs = [ - ("ACGTACGTACGTACGT", "ACGTACGTACGTACGT"), - ("AAAAAAAAAAAAAAAA", "CCCCCCCCCCCCCCCC"), - ("ATCGATCGATCGATCG", "TAGCTAGCTAGCTAGC"), - ("ACGTACGTACGTACGTACGTACGT", "TTTTAAAACCCCGGGG"), - ("AACCGGTTAACCGGTT", "GGTTAACCGGTTAACC"), - ]; - - let ctx = DistanceContext::new(4); - let cpp = crate::original_ffi::OriginalSpearmanDist::new(4); - - for (seq1, seq2) in &test_pairs { - let rust_d = spearman_dist_seq(&ctx, seq1.as_bytes(), seq2.as_bytes()); - let cpp_d = cpp.get_dist_seq(seq1, seq2); - - let diff = (rust_d - cpp_d).abs(); - assert!( - diff < 1e-10, - "SpearmanDist seq mismatch for '{seq1}' vs '{seq2}': rust={rust_d} cpp={cpp_d}" - ); - } - } - - #[test] - fn ffi_spearman_dist_profile_equivalence() { - let ctx = DistanceContext::new(4); - let cpp = crate::original_ffi::OriginalSpearmanDist::new(4); - - let rust_kmap = KmerMap::new(4, true); - let seqs = ["ACGTACGTACGTACGT", "AAAAAAAAAAAAAAAA", "ATCGATCGATCGATCG"]; - - let profiles: Vec> = seqs - .iter() - .map(|s| { - Profiler::new(4, s.as_bytes(), &rust_kmap) - .get_profile() - .to_vec() - }) - .collect(); - - for i in 0..profiles.len() { - for j in (i + 1)..profiles.len() { - let rust_d = spearman_dist_profiles(&ctx, &profiles[i], &profiles[j]); - let cpp_d = cpp.get_dist_profile(&profiles[i], &profiles[j]); - - let diff = (rust_d - cpp_d).abs(); - assert!( - diff < 1e-10, - "SpearmanDist profile mismatch: rust={rust_d} cpp={cpp_d}" - ); - } - } - } - - #[test] - fn ffi_spearman_normalized_equivalence() { - let mut ctx = DistanceContext::new(4); - ctx.set_normalization(true); - let cpp = crate::original_ffi::OriginalSpearmanDist::new(4); - cpp.set_normalization(true); - - let seq1 = "ACGTACGTACGTACGT"; - let seq2 = "AAAAAAAAAAAAAAAA"; - - let rust_d = spearman_dist_seq(&ctx, seq1.as_bytes(), seq2.as_bytes()); - let cpp_d = cpp.get_dist_seq(seq1, seq2); - - let diff = (rust_d - cpp_d).abs(); - assert!( - diff < 1e-10, - "SpearmanDist normalized mismatch: rust={rust_d} cpp={cpp_d}" - ); - } } diff --git a/crates/maxbin-rs/src/emanager.rs b/crates/maxbin-rs/src/emanager.rs index d33f68c..f6d1c32 100644 --- a/crates/maxbin-rs/src/emanager.rs +++ b/crates/maxbin-rs/src/emanager.rs @@ -774,7 +774,7 @@ mod tests { // NOTE: get_prob_abund_matches_cpp was removed — it computed the "expected" // value using the same lgamma and formula as the function under test, making - // it tautological. Real coverage comes from proptest_emanager.rs (algebraic - // properties on the real function) and emanager_equivalence.rs (full pipeline - // byte-for-byte comparison against C++). + // it tautological. Real coverage lives in the maxbin-rs-equivalence crate: + // its proptest_emanager.rs (algebraic properties on the real function) and + // emanager_equivalence.rs (full pipeline byte-for-byte comparison against C++). } diff --git a/crates/maxbin-rs/src/fasta.rs b/crates/maxbin-rs/src/fasta.rs index f3ac041..46fb612 100644 --- a/crates/maxbin-rs/src/fasta.rs +++ b/crates/maxbin-rs/src/fasta.rs @@ -130,69 +130,6 @@ mod tests { assert_eq!(records[0].seq, b"ACGTTGCA"); } - /// Decompress a .gz file to a temp path for the C++ parser (which can't read gzip). - fn decompress_gz(gz_path: &std::path::Path) -> std::path::PathBuf { - use std::io::{Read, Write}; - let out_path = std::env::temp_dir().join("maxbin_rs_test_contigs.fa"); - let file = std::fs::File::open(gz_path).unwrap(); - let mut decoder = flate2::read::GzDecoder::new(file); - let mut contents = Vec::new(); - decoder.read_to_end(&mut contents).unwrap(); - let mut out = std::fs::File::create(&out_path).unwrap(); - out.write_all(&contents).unwrap(); - out_path - } - - #[test] - fn equivalence_with_original_fasta_parser() { - // This test requires nf-core test data (only available in devshell). - let gz_path = match std::env::var("MAXBIN2_TEST_CONTIGS") { - Ok(p) => std::path::PathBuf::from(p), - Err(_) => return, - }; - - // Parse with our Rust implementation (handles gzip natively) - let rust_records = parse_file(&gz_path).expect("rust parser failed"); - - // Decompress for the C++ parser, then parse with the original - let plain_path = decompress_gz(&gz_path); - let original = crate::original_ffi::OriginalFastaReader::new(&plain_path); - - // Same number of records - assert_eq!( - rust_records.len(), - original.num_records() as usize, - "record count mismatch" - ); - - // Field-by-field comparison for every record - for (i, rust_rec) in rust_records.iter().enumerate() { - let idx = i as u32; - - assert_eq!( - rust_rec.header, - original.header(idx), - "header mismatch at record {i}" - ); - - assert_eq!( - rust_rec.seq, - original.seq(idx), - "sequence mismatch at record {i} (header: {})", - rust_rec.header - ); - - assert_eq!( - rust_rec.len(), - original.seq_len(idx) as usize, - "length mismatch at record {i} (header: {})", - rust_rec.header - ); - } - - std::fs::remove_file(&plain_path).ok(); - } - #[test] fn parse_nfcore_test_data() { // This test uses the nf-core test contigs if available via env var. diff --git a/crates/maxbin-rs/src/kmer_map.rs b/crates/maxbin-rs/src/kmer_map.rs index d69a972..ae1aedb 100644 --- a/crates/maxbin-rs/src/kmer_map.rs +++ b/crates/maxbin-rs/src/kmer_map.rs @@ -249,53 +249,4 @@ mod tests { let rev = km.get_reverse_mapping_idx(0); assert_eq!(rev, km.get_mapping(b"TTTT")); } - - #[test] - fn ffi_equivalence_kmer4_symmetric() { - let rust_km = KmerMap::new(4, true); - let cpp_km = crate::original_ffi::OriginalKmerMap::new(4, true); - - assert_eq!(rust_km.get_entry_num() as i32, cpp_km.entry_num()); - - let cpp_table = cpp_km.get_full_table(); - let rust_table = rust_km.get_full_table(); - - assert_eq!(rust_table.len(), cpp_table.len()); - for i in 0..rust_table.len() { - assert_eq!(rust_table[i], cpp_table[i], "mapping mismatch at index {i}"); - } - } - - #[test] - fn ffi_equivalence_kmer4_nonsymmetric() { - let rust_km = KmerMap::new(4, false); - let cpp_km = crate::original_ffi::OriginalKmerMap::new(4, false); - - assert_eq!(rust_km.get_entry_num() as i32, cpp_km.entry_num()); - - let cpp_table = cpp_km.get_full_table(); - let rust_table = rust_km.get_full_table(); - - for i in 0..rust_table.len() { - assert_eq!(rust_table[i], cpp_table[i]); - } - } - - #[test] - fn ffi_equivalence_kmer3_symmetric() { - let rust_km = KmerMap::new(3, true); - let cpp_km = crate::original_ffi::OriginalKmerMap::new(3, true); - - assert_eq!(rust_km.get_entry_num() as i32, cpp_km.entry_num()); - - let cpp_table = cpp_km.get_full_table(); - let rust_table = rust_km.get_full_table(); - - for i in 0..rust_table.len() { - assert_eq!( - rust_table[i], cpp_table[i], - "kmer3 symmetric mapping mismatch at index {i}" - ); - } - } } diff --git a/crates/maxbin-rs/src/lib.rs b/crates/maxbin-rs/src/lib.rs index d37c209..7ac1606 100644 --- a/crates/maxbin-rs/src/lib.rs +++ b/crates/maxbin-rs/src/lib.rs @@ -10,7 +10,3 @@ pub mod paths; pub mod pipeline; pub mod profiler; pub mod quicksort; - -// Always compiled — the C++ FFI is linked by build.rs regardless, and -// integration tests need access to these bindings. -pub mod original_ffi; diff --git a/crates/maxbin-rs/src/main.rs b/crates/maxbin-rs/src/main.rs index 351ac19..e6c8a55 100644 --- a/crates/maxbin-rs/src/main.rs +++ b/crates/maxbin-rs/src/main.rs @@ -5,7 +5,6 @@ fn main() { maxbin_rs::cli::Command::Filter(ref args) => maxbin_rs::pipeline::run_filter(args), maxbin_rs::cli::Command::Seeds(ref args) => maxbin_rs::pipeline::run_seeds(args), maxbin_rs::cli::Command::Em(ref args) => maxbin_rs::pipeline::run_em(args), - maxbin_rs::cli::Command::CppEm(ref args) => maxbin_rs::pipeline::run_cpp_em(args), maxbin_rs::cli::Command::SamToAbund(ref args) => { maxbin_rs::pipeline::run_sam_to_abund(args) } diff --git a/crates/maxbin-rs/src/normal_distribution.rs b/crates/maxbin-rs/src/normal_distribution.rs index 90106bb..88abc3e 100644 --- a/crates/maxbin-rs/src/normal_distribution.rs +++ b/crates/maxbin-rs/src/normal_distribution.rs @@ -66,37 +66,4 @@ mod tests { // They should differ slightly assert!(our_val != std_pi_val); } - - #[test] - fn ffi_equivalence() { - let test_cases = [ - (0.0, 1.0, 0.0), - (0.0, 1.0, 1.0), - (0.0, 1.0, -1.0), - (5.0, 2.0, 5.0), - (5.0, 2.0, 3.0), - (5.0, 2.0, 7.0), - (10.0, 0.5, 10.0), - (10.0, 0.5, 9.5), - (-3.0, 1.5, -3.0), - (-3.0, 1.5, 0.0), - (0.0, 1.0, 3.5), - (100.0, 10.0, 95.0), - ]; - - for &(mean, std, input) in &test_cases { - let rust_nd = NormalDistribution::new(mean, std); - let cpp_nd = crate::original_ffi::OriginalNormalDistribution::new(mean, std); - - let rust_p = rust_nd.prob(input); - let cpp_p = cpp_nd.prob(input); - - assert_eq!( - rust_p.to_bits(), - cpp_p.to_bits(), - "NormalDist not bit-identical for mean={mean}, std={std}, input={input}: \ - rust={rust_p:e} cpp={cpp_p:e}" - ); - } - } } diff --git a/crates/maxbin-rs/src/pipeline.rs b/crates/maxbin-rs/src/pipeline.rs index e019205..ea22eee 100644 --- a/crates/maxbin-rs/src/pipeline.rs +++ b/crates/maxbin-rs/src/pipeline.rs @@ -8,7 +8,7 @@ /// 5. Recursive binning: seed → EM → re-seed output bins → repeat (up to 5 levels) /// 6. Merge noclass, filter small bins, sort by abundance /// 7. Write output files -use crate::cli::{CppEmArgs, EmArgs, FilterArgs, PipelineArgs, SamToAbundArgs, SeedsArgs}; +use crate::cli::{EmArgs, FilterArgs, PipelineArgs, SamToAbundArgs, SeedsArgs}; use crate::fasta; use std::io::Write; use std::path::{Path, PathBuf}; @@ -348,26 +348,6 @@ pub fn run_em(args: &EmArgs) -> Result<(), String> { Ok(()) } -/// Run the `cpp-em` subcommand: run the original C++ EM via FFI. -/// The C++ EManager handles its own data loading, EM, classification, -/// and file writing internally — we time the whole thing as one block. -pub fn run_cpp_em(args: &CppEmArgs) -> Result<(), String> { - eprintln!("Running original C++ EM via FFI..."); - eprintln!(" contigs: {}", args.contig.display()); - eprintln!(" abund: {}", args.abund.display()); - eprintln!(" seed: {}", args.seed.display()); - eprintln!(" threads: {}", args.thread); - - let t = std::time::Instant::now(); - let em = crate::original_ffi::OriginalEManager::new(&args.contig, &args.abund, &args.out); - em.set_thread_num(args.thread as i32); - let result = em.run(&args.seed); - let secs = t.elapsed().as_secs_f64(); - - eprintln!("C++ EM completed in {secs:.1}s (returned: {result})"); - Ok(()) -} - /// Run the `sam-to-abund` subcommand: compute abundance from a SAM file. pub fn run_sam_to_abund(args: &SamToAbundArgs) -> Result<(), String> { compute_abundance_from_sam(&args.sam, &args.out)?; diff --git a/crates/maxbin-rs/src/profiler.rs b/crates/maxbin-rs/src/profiler.rs index 2b92937..f3afafa 100644 --- a/crates/maxbin-rs/src/profiler.rs +++ b/crates/maxbin-rs/src/profiler.rs @@ -148,45 +148,4 @@ mod tests { let prof = Profiler::new(4, seq, &kmap); assert!(prof.get_profile().iter().all(|&v| v == 0.0)); } - - #[test] - fn ffi_equivalence() { - let test_seqs = [ - "ACGTACGTACGTACGT", - "AAAAAAAAAAAA", - "ACGTTTTTGGGGCCCC", - "ATCGATCGATCGATCGATCG", - "ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGT", - ]; - - let rust_kmap = KmerMap::new(4, true); - let cpp_kmap = crate::original_ffi::OriginalKmerMap::new(4, true); - let entry_num = rust_kmap.get_entry_num(); - - for seq in &test_seqs { - let rust_prof = Profiler::new(4, seq.as_bytes(), &rust_kmap); - let cpp_prof = crate::original_ffi::OriginalProfiler::new(4, seq, &cpp_kmap); - - let rust_profile = rust_prof.get_profile(); - let cpp_profile = cpp_prof.get_profile(entry_num as i32); - - assert_eq!(rust_profile.len(), cpp_profile.len()); - for i in 0..rust_profile.len() { - assert_eq!( - rust_profile[i].to_bits(), - cpp_profile[i].to_bits(), - "profile not bit-identical at index {i} for seq '{seq}': rust={:e} cpp={:e}", - rust_profile[i], - cpp_profile[i] - ); - } - - let rust_pn = rust_prof.get_percent_n(); - let cpp_pn = cpp_prof.get_percent_n(); - assert!( - (rust_pn - cpp_pn).abs() < 1e-6, - "percent_n mismatch for seq '{seq}': rust={rust_pn} cpp={cpp_pn}" - ); - } - } } diff --git a/crates/maxbin-rs/src/quicksort.rs b/crates/maxbin-rs/src/quicksort.rs index c1b4530..20ca96a 100644 --- a/crates/maxbin-rs/src/quicksort.rs +++ b/crates/maxbin-rs/src/quicksort.rs @@ -111,41 +111,4 @@ mod tests { sort_descending(&mut arr, None); assert_eq!(arr, vec![3.0, 3.0, 2.0, 1.0, 1.0]); } - - #[test] - fn equivalence_with_original() { - // Test a variety of inputs against the C++ quickSort via FFI. - let test_cases: Vec> = vec![ - vec![1.0, 5.0, 3.0, 2.0, 4.0], - vec![5.0, 4.0, 3.0, 2.0, 1.0], - vec![1.0, 1.0, 1.0], - vec![3.0, 1.0, 4.0, 1.0, 5.0, 9.0, 2.0, 6.0, 5.0, 3.0, 5.0], - vec![42.0], - vec![1.0, 2.0], - vec![2.0, 1.0], - vec![0.001, 0.01, 0.1, 1.0, 10.0, 100.0], - vec![-1.0, 0.0, 1.0, -0.5, 0.5], - ]; - - for (case_idx, case) in test_cases.iter().enumerate() { - // Rust version - let mut rust_arr = case.clone(); - let mut rust_idx: Vec = (0..case.len() as i32).collect(); - sort_descending(&mut rust_arr, Some(&mut rust_idx)); - - // Original C++ version - let mut orig_arr = case.clone(); - let mut orig_idx: Vec = (0..case.len() as i32).collect(); - crate::original_ffi::original_quicksort(&mut orig_arr, Some(&mut orig_idx)); - - assert_eq!( - rust_arr, orig_arr, - "array mismatch for case {case_idx}: {case:?}" - ); - assert_eq!( - rust_idx, orig_idx, - "index mismatch for case {case_idx}: {case:?}" - ); - } - } } diff --git a/crates/maxbin-rs/tests/cli-integration.sh b/crates/maxbin-rs/tests/cli-integration.sh index d6e7acf..2c7c293 100644 --- a/crates/maxbin-rs/tests/cli-integration.sh +++ b/crates/maxbin-rs/tests/cli-integration.sh @@ -301,26 +301,6 @@ assert_exit_code "4.9 em-missing-out" 1 \ echo "" -# ========================================================================= -# 5. cpp-em subcommand -# ========================================================================= -echo "--- 5. cpp-em subcommand ---" - -D="$WORK/cpp_em" - -# 5.1 basic -mkdir -p "$D/basic" -assert_exit_code "5.1 cpp-em-basic" 0 \ - maxbin-rs cpp-em --contig "$FILTERED" --abund "$ABUND" --seed "$SEED" --out "$D/basic/test" -assert_file_exists "$D/basic/test.summary" "5.1 summary" - -# 5.2 underscore subcommand variant -mkdir -p "$D/underscore" -assert_exit_code "5.2 cpp_em-underscore" 0 \ - maxbin-rs cpp_em --contig "$FILTERED" --abund "$ABUND" --seed "$SEED" --out "$D/underscore/test" - -echo "" - # ========================================================================= # 6. sam-to-abund subcommand # ========================================================================= diff --git a/crates/maxbin-rs/tests/pipeline-stages.sh b/crates/maxbin-rs/tests/pipeline-stages.sh index 7588ebb..2292cba 100755 --- a/crates/maxbin-rs/tests/pipeline-stages.sh +++ b/crates/maxbin-rs/tests/pipeline-stages.sh @@ -308,10 +308,13 @@ STAGE4_START=$SECONDS CPP_DIR=$(mktemp -d) RUST_DIR=$(mktemp -d) -# Subcommands: run C++ EM and Rust EM separately. Each reports its own timing. +# Run C++ EM and Rust EM separately. Each reports its own timing. +# The C++ EM is the maxbin-rs-cpp-em binary from the maxbin-rs-equivalence +# crate (it links the original MaxBin2 C++ via FFI); it used to be the +# `maxbin-rs cpp-em` subcommand before the FFI was split into its own crate. # Both use prob_threshold=0.5 to match the C++ FFI's hardcoded default. # (Our default changed to 0.9 in v0.2, but the FFI can't be overridden.) -maxbin-rs cpp-em --contig "$FILTERED_CONTIGS" --abund "$ORIG_ABUND" \ +maxbin-rs-cpp-em --contig "$FILTERED_CONTIGS" --abund "$ORIG_ABUND" \ --seed "$ORIG_SEED" --out "$CPP_DIR/test" --thread 1 maxbin-rs em --contig "$FILTERED_CONTIGS" --abund "$ORIG_ABUND" \ diff --git a/flake.nix b/flake.nix index 9124d19..88e8608 100644 --- a/flake.nix +++ b/flake.nix @@ -75,6 +75,11 @@ path: type: (craneLib.filterCargoSources path type) || builtins.match ".*/(nix|vendor|tests)/.*" path != null; }; + # Shared build arguments. MAXBIN2_SRC_TARBALL is consumed by the + # maxbin-rs-equivalence crate's build.rs (it extracts and links the + # original MaxBin2 C++ for FFI equivalence testing). The core + # maxbin-rs crate has no build.rs and ignores it; the workspace-wide + # clippy and nextest checks need it because they span both crates. commonArgs = { inherit src; pname = "maxbin-rs"; @@ -84,6 +89,13 @@ }; cargoArtifacts = craneLib.buildDepsOnly commonArgs; + # Build arguments for the core crate alone — no C++ toolchain or + # tarball needed. The published binaries (maxbin-rs, maxbin-rs-em) + # build only this crate so they carry zero C/C++. + coreArgs = builtins.removeAttrs commonArgs [ "MAXBIN2_SRC_TARBALL" ] // { + cargoExtraArgs = "-p maxbin-rs"; + }; + # ===================================================================== # Packages # ===================================================================== @@ -137,7 +149,7 @@ # Supports `em`, `filter`, and `seeds` subcommands only. # Pipeline mode errors out telling the user to use the full package. maxbin-rs-em = craneLib.buildPackage ( - commonArgs + coreArgs // { inherit cargoArtifacts; pname = "maxbin-rs-em"; @@ -146,7 +158,7 @@ # The Rust reimplementation — this is the main output of this project. maxbin-rs = craneLib.buildPackage ( - commonArgs + coreArgs // { inherit cargoArtifacts; # Wrap the binary so HMMER, Bowtie2, and FragGeneScan are on $PATH @@ -172,23 +184,27 @@ } ); - # Same as maxbin-rs but with C++ LTO enabled for the FFI library. - # Used for A/B benchmarking to test whether -flto closes the ~8x gap. - maxbin-rs-lto = craneLib.buildPackage ( + # The equivalence crate's cpp-em binary with C++ LTO enabled for the + # FFI library. Used for A/B benchmarking to test whether -flto closes + # the performance gap. Produces $out/bin/maxbin-rs-cpp-em. + maxbin-rs-cpp-em-lto = craneLib.buildPackage ( commonArgs // { inherit cargoArtifacts; - pname = "maxbin-rs-lto"; + pname = "maxbin-rs-cpp-em-lto"; + cargoExtraArgs = "-p maxbin-rs-equivalence --bin maxbin-rs-cpp-em"; MAXBIN2_CPP_LTO = "1"; - postInstall = '' - wrapProgram $out/bin/maxbin-rs \ - --prefix PATH : "${ - pkgs.lib.makeBinPath [ - pkgs.hmmer - pkgs.bowtie2 - ] - }:${fraggenescan}/libexec/FragGeneScan" - ''; + } + ); + + # The equivalence crate's cpp-em binary (no LTO), the baseline for the + # A/B benchmark above. Produces $out/bin/maxbin-rs-cpp-em. + maxbin-rs-cpp-em = craneLib.buildPackage ( + commonArgs + // { + inherit cargoArtifacts; + pname = "maxbin-rs-cpp-em"; + cargoExtraArgs = "-p maxbin-rs-equivalence --bin maxbin-rs-cpp-em"; } ); @@ -298,12 +314,15 @@ maxbin2-trace maxbin2-f64-trace maxbin-rs - maxbin-rs-lto + maxbin-rs-em + maxbin-rs-cpp-em + maxbin-rs-cpp-em-lto fraggenescan fraggenescan-rs rust ; inherit (pkgs) cargo-nextest; + inherit maxbin2-src-tarball; inherit datasets intermediates; }; diff --git a/nix/tests.nix b/nix/tests.nix index 8c1cd97..6372210 100644 --- a/nix/tests.nix +++ b/nix/tests.nix @@ -20,11 +20,14 @@ maxbin2-trace, maxbin2-f64-trace, maxbin-rs, - maxbin-rs-lto, + maxbin-rs-em, + maxbin-rs-cpp-em, + maxbin-rs-cpp-em-lto, fraggenescan, fraggenescan-rs, rust, cargo-nextest, + maxbin2-src-tarball, datasets, intermediates, }: @@ -175,7 +178,7 @@ in ]; } '' - export MAXBIN2_SRC_TARBALL="${maxbin-rs.MAXBIN2_SRC_TARBALL}" + export MAXBIN2_SRC_TARBALL="${maxbin2-src-tarball}" cp -r ${../.}/* . chmod -R u+w . cargo nextest run emanager_precision_divergence --no-fail-fast @@ -190,7 +193,8 @@ in runtimeInputs = [ maxbin2 maxbin-rs - maxbin-rs-lto + maxbin-rs-cpp-em + maxbin-rs-cpp-em-lto ]; text = '' echo "=== C++ LTO Benchmark: CAMI I High EM ===" @@ -203,7 +207,7 @@ in trap 'rm -rf "$WORK"' EXIT echo "--- Filtering contigs ---" - "${maxbin-rs}/bin/maxbin-rs" filter -contig "$CONTIGS" -out "$WORK/test" + "${maxbin-rs}/bin/maxbin-rs" filter --contig "$CONTIGS" --out "$WORK/test" FILTERED="$WORK/test.contig.tmp" ABUND="$INTERMEDIATES/abund" @@ -212,25 +216,25 @@ in echo "" echo "--- C++ EM (baseline, no LTO) ---" CPP_DIR=$(mktemp -d) - MAXBIN_RS_DETERMINISTIC=1 "${maxbin-rs}/bin/maxbin-rs" cpp-em \ - -contig "$FILTERED" -abund "$ABUND" -seed "$SEED" \ - -out "$CPP_DIR/test" -thread 1 + MAXBIN_RS_DETERMINISTIC=1 "${maxbin-rs-cpp-em}/bin/maxbin-rs-cpp-em" \ + --contig "$FILTERED" --abund "$ABUND" --seed "$SEED" \ + --out "$CPP_DIR/test" --thread 1 rm -rf "$CPP_DIR" echo "" echo "--- C++ EM (with -flto) ---" LTO_DIR=$(mktemp -d) - MAXBIN_RS_DETERMINISTIC=1 "${maxbin-rs-lto}/bin/maxbin-rs" cpp-em \ - -contig "$FILTERED" -abund "$ABUND" -seed "$SEED" \ - -out "$LTO_DIR/test" -thread 1 + MAXBIN_RS_DETERMINISTIC=1 "${maxbin-rs-cpp-em-lto}/bin/maxbin-rs-cpp-em" \ + --contig "$FILTERED" --abund "$ABUND" --seed "$SEED" \ + --out "$LTO_DIR/test" --thread 1 rm -rf "$LTO_DIR" echo "" echo "--- Rust EM ---" RUST_DIR=$(mktemp -d) MAXBIN_RS_DETERMINISTIC=1 "${maxbin-rs}/bin/maxbin-rs" em \ - -contig "$FILTERED" -abund "$ABUND" -seed "$SEED" \ - -out "$RUST_DIR/test" -thread 1 + --contig "$FILTERED" --abund "$ABUND" --seed "$SEED" \ + --out "$RUST_DIR/test" --thread 1 rm -rf "$RUST_DIR" echo "" @@ -249,6 +253,7 @@ in ]; text = '' export MAXBIN2_TEST_CONTIGS="${datasets.bfragilis.contigs}" + export MAXBIN2_SRC_TARBALL="${maxbin2-src-tarball}" # Copy source to a writable temp dir (Nix store is read-only). WORK=$(mktemp -d) cp -r ${../.}/* "$WORK/" @@ -340,28 +345,33 @@ in echo " $label: $total lines, $calls calls, $vector vector, $scalar scalar" } - # disasm_variant BINARY OUTDIR VARIANT_LABEL - # Disassemble all four EM hot functions from BINARY into OUTDIR. + # disasm_variant CPP_BINARY RUST_BINARY OUTDIR VARIANT_LABEL + # Disassemble the four EM hot functions into OUTDIR. The C++ functions + # come from CPP_BINARY (the cpp-em binary, which is what the standard + # vs LTO A/B varies); the Rust functions come from RUST_BINARY. They + # have to be pulled from different binaries: the cpp-em binary never + # calls the Rust EM, so those symbols are dead-code-eliminated from it, + # while the maxbin-rs binary carries no C++ EManager at all. disasm_variant() { - local binary="$1" outdir="$2" label="$3" + local cpp_binary="$1" rust_binary="$2" outdir="$3" label="$4" mkdir -p "$outdir" - disasm_function "$binary" \ + disasm_function "$cpp_binary" \ 'T EManager::run_EM\(int\)$' \ "$outdir/cpp-run_EM.s" \ "$label — C++ EManager::run_EM(int)" - disasm_function "$binary" \ + disasm_function "$cpp_binary" \ 'EManager::get_prob_abund\(double,? ?double\)' \ "$outdir/cpp-get_prob_abund.s" \ "$label — C++ EManager::get_prob_abund(double, double)" - disasm_function "$binary" \ + disasm_function "$rust_binary" \ 'T maxbin_rs.*emanager.*run_em$' \ "$outdir/rust-run_em.s" \ "$label — Rust run_em" - disasm_function "$binary" \ + disasm_function "$rust_binary" \ 't maxbin_rs.*emanager.*compute_abund_prob_for_contig$' \ "$outdir/rust-compute_abund.s" \ "$label — Rust compute_abund_prob_for_contig" @@ -379,11 +389,6 @@ in } ''; - # Resolve the unwrapped binary (skip the shell wrapper). - binaryOf = pkg: "${pkg}/bin/.maxbin-rs-wrapped"; - # Fallback if the wrapper layout differs. - binaryFallback = pkg: "${pkg}/bin/maxbin-rs"; - in runCommand "disasm-em" { @@ -392,26 +397,27 @@ in '' ${disasmLib} - # Find the actual ELF binary (wrapProgram creates a shell wrapper). - find_binary() { - local pkg="$1" - if [ -f "$pkg/bin/.maxbin-rs-wrapped" ]; then - echo "$pkg/bin/.maxbin-rs-wrapped" - else - echo "$pkg/bin/maxbin-rs" - fi - } - - STD_BIN=$(find_binary "${maxbin-rs}") - LTO_BIN=$(find_binary "${maxbin-rs-lto}") + # The C++ EManager symbols come from the equivalence crate's cpp-em + # binary (which links the original MaxBin2 C++ via FFI); the standard + # vs LTO A/B varies only this binary (MAXBIN2_CPP_LTO=1). The Rust EM + # symbols come from the maxbin-rs binary instead: the cpp-em binary + # only calls the C++ EM, so the Rust EM is dead-code-eliminated from + # it. The Rust codegen is identical across the two columns (LTO affects + # only the C++ FFI), so the same RUST_BIN feeds both. + STD_BIN="${maxbin-rs-cpp-em}/bin/maxbin-rs-cpp-em" + LTO_BIN="${maxbin-rs-cpp-em-lto}/bin/maxbin-rs-cpp-em" + # maxbin-rs-em is the unwrapped full binary (no wrapProgram shell + # wrapper), and its pipeline actually calls the Rust EM, so run_em and + # compute_abund_prob_for_contig survive in it. + RUST_BIN="${maxbin-rs-em}/bin/maxbin-rs" mkdir -p $out/standard $out/unity - echo "Disassembling standard build: $STD_BIN" - disasm_variant "$STD_BIN" "$out/standard" "standard" + echo "Disassembling standard build: $STD_BIN (Rust from $RUST_BIN)" + disasm_variant "$STD_BIN" "$RUST_BIN" "$out/standard" "standard" - echo "Disassembling unity/LTO build: $LTO_BIN" - disasm_variant "$LTO_BIN" "$out/unity" "unity (MAXBIN2_CPP_LTO=1)" + echo "Disassembling unity/LTO build: $LTO_BIN (Rust from $RUST_BIN)" + disasm_variant "$LTO_BIN" "$RUST_BIN" "$out/unity" "unity (MAXBIN2_CPP_LTO=1)" # Side-by-side comparison {