HYDRA v2.0 (hydra_v2_sa)

Software Overview

Maximum Parsimony Phylogenetic Reconstruction

Given an aligned matrix of discrete morphological or molecular characters across n operational taxonomic units (taxa), hydra_v2_sa searches for a rooted binary phylogenetic tree topology that minimizes the total evolutionary mutation count (parsimony length) required to explain the observed leaf sequences.

1. Alignment Compression Filters invariant & singleton sites into an exact base_cost; sorts informative columns by decreasing variability.

2. Adaptive Initialization Uniform random trees for n < 200; AVX2-accelerated randomized greedy 1stRotuGbr stepwise addition for n ≥ 200.

3. Compound-SPR Annealing Goloboff two-pass state buffers + store-free AVX2 kernels with exact Metropolis early-exit bounds (max_delta).

4. Elite Reheating & Polish Dynamic elite reheating (β = 1.4) on stagnation, followed by exhaustive strict SPR local-search descent and exact verification.

Architecture & Algorithms

Core Software Characteristics

hydra_v2_sa pairs domain-specific phylogenetic tree invariants with low-level SIMD micro-optimizations to evaluate millions of Subtree Pruning and Regrafting (SPR) rearrangements per second.

01 • DATA ENCODING

Bitmask Encoding & Variability Sorting

Each character state is encoded as an 8-bit bitmask supporting standard morphological states (0–3), nucleotides (A, C, G, T/U), gaps as a 5th state (-), missing data (?, N, X), full IUPAC ambiguity codes (R, Y, M, K, S, W, B, D, H, V), and bracketed polymorphic states (e.g., [01], [AC]).

Zero-Loss Compression: Invariant sites (0 steps) and parsimony-uninformative singletons are stripped upfront; singleton steps are accumulated into an exact constant base_cost.
Early-Exit Ordering: Informative sites are sorted by decreasing state diversity so mutations cluster in the leading 32-byte SIMD blocks.

02 • GOLOBOFF + AVX2 ENGINE

Two-Pass PhyloTree & Store-Free SIMD Kernels

Internal and leaf sequences are stored in 32-byte aligned buffers (S_down, S_up, S_edge) implementing Pablo Goloboff's two-pass parsimony framework.

Selective Top-Down Pass: When pruning subtree p, HYDRA flags eligible target edges (need_edge, need_up) and computes top-down states only where needed on the residual tree T_R.
O(1) Branch Reinsertion: Evaluates 32 sites per CPU cycle in registers via 256-bit AVX2 bitwise masks and hardware POPCNT, aborting immediately whenever accumulated steps exceed max_delta.

03 • SIMULATED ANNEALING

Compound SPR Neighborhood & Metropolis Pruning

At temperature T, for a candidate pruned subtree p (drawn without immediate repetition in shuffled epochs), HYDRA scans all valid regrafting branches v ≠ s on T_R to select the best compound neighbor.

Exact Metropolis Bound: Drawing u ∈ (0, 1), a cost change Δf satisfies e^−Δf/T > u iff Δf < −T ln(u), yielding an exact integer ceiling max_degrad = ⌈-T ln(u)⌉ - 1 passed directly to the AVX2 kernel.
Dynamic Reheating & Intensification: Stagnation over 18 plateaus triggers elite reheating (T ← max(1.4T, 0.35T_i·0.85^r-1)), ending with strict SPR local-search descent (max_degrad = -1).

04 • EXECUTION MODES

Sequential Default & Cooperative OpenMP Chains

HYDRA adapts seamlessly from lightweight single-core scripting to multi-core workstation runs while automatically tuning hyperparameters around the n = 200 taxa threshold.

Sequential Mode (Default): Single deterministic Markov chain delivering sub-second to ~1.3s runtimes on medium alignments and ~4.4s on the 500-taxon zilla.txt dataset.
Parallel Mode (--parallel): Spawns independent OpenMP Markov chains from distinct seeds/topologies that synchronize and share the global elite tree at plateau boundaries and reheating phases.

Complete CLI Reference

Command-Line Options & Default Values

Invocation syntax: ./hydra [options] <instance_file> (where <instance_file> may be passed positionally or via -i / --input). Automatic defaults adapt dynamically to the number of taxa n (with total tree nodes N = 2n − 1).

Execution Mode & General Options

10 Options

Option Flag	Argument	Default Value	Description & Automatic Resolution Behavior
--parallel	—	false (disabled)	Enables multi-chain Parallel Simulated Annealing using OpenMP threads. By default, HYDRA runs a single non-parallel (sequential) Simulated Annealing Markov chain.
-j, --threads	<int>	0 (auto)	Number of concurrent OpenMP Markov chains when `--parallel` is enabled. When set to `0`, automatically resolves to `4 × hardware_threads` if n < 200 taxa, or `1 × hardware_threads` if n ≥ 200 taxa (`1` in sequential mode).
-i, --input	<file>	Positional arg	Path to the input character alignment matrix in PHYLIP-like (`.txt`) or TNT (`.sg`) format. Can also be supplied directly as a positional argument.
-o, --output	<file>	"" (none)	Destination file path where the best discovered phylogenetic tree topology is exported in standard Newick format.
-t, --time-limit	<float>	0 (unlimited)	Maximum wall-clock execution time in seconds. When `≤ 0`, execution runs until the cooling/reheating schedule and final local-search intensification complete.
-s, --seed	<uint>	8	Base seed for the 32-bit Mersenne Twister PRNG (`std::mt19937`). Under `--parallel`, each chain `tid` is deterministically seeded with `seed + tid × 13337`.
-r, --random-init	—	auto (if n < 200)	Forces uniform random initial tree topologies regardless of instance size, maximizing multi-basin structural diversity across runs and Markov chains.
-g, --greedy-init	—	auto (if n >= 200)	Forces randomized greedy stepwise taxon addition (`1stRotuGbr`, randomized tie/candidate factor `0.90`) accelerated to O(n²·k) via two-pass AVX2 evaluation.
-q, --quiet	—	false (verbose)	Quiet mode: suppresses configuration banners, plateau logs, and reheating traces, printing solely the verified final integer parsimony score to `stdout`.
-h, --help	—	—	Displays the command-line usage summary and option reference, then exits immediately with status code `0`.

Simulated Annealing Schedule Parameters

5 Options

Option Flag	Argument	Default Value	Description & Automatic Resolution Behavior
-T, --temp-init	<float>	6.0	Initial Simulated Annealing temperature T_i (clamped internally to `≥ 10^-4`). Controls initial Metropolis acceptance probability for uphill SPR moves.
-F, --temp-final	<float>	0.005	Final stopping temperature T_f (clamped internally to `≥ 10^-6`). Cooling terminates once T < T_f and all dynamic reheats are exhausted.
-c, --cooling	<float>	0.91	Geometric cooling factor γ ∈ (0, 1) applied at the end of each temperature plateau (T ← γ · T). Out-of-range values fall back to `0.91`.
-L, --mcl	<int>	0 (auto)	Markov chain plateau length (compound SPR trials per temperature step). When `0`, resolves automatically to `2000` if n < 200; for n ≥ 200, resolves to `max(600, ⌊3(2n-1)/5⌋)` with `--parallel` or `max(800, ⌊4(2n-1)/5⌋)` sequentially.
-R, --reheats	<int>	auto (20 \| 12)	Maximum number of dynamic reheating cycles triggered after 18 stagnant plateaus or upon reaching T_f. Automatically set to `20` if n < 200 or `12` if n ≥ 200 (set `0` to disable reheating).

Automatic Regime for Small / Medium Instances (n < 200)

Initial Topology: Uniform random binary tree (-r) to maximize basin exploration.
Markov Plateau Length (-L): 2,000 compound SPR iterations per plateau.
Max Reheats (-R): 20 elite reheating cycles.
OpenMP Chains (--parallel -j 0): 4 × hardware_threads concurrent chains.

Automatic Regime for Large Instances (n ≥ 200, e.g. `zilla.txt`)

Initial Topology: Randomized greedy stepwise addition 1stRotuGbr (-g).
Markov Plateau Length (-L): max(800, ⌊0.8(2n-1)⌋) sequentially; max(600, ⌊0.6(2n-1)⌋) with --parallel.
Max Reheats (-R): 12 elite reheating cycles, with early low-temperature stagnation exit (T ≤ 0.05).
OpenMP Chains (--parallel -j 0): 1 × hardware_threads concurrent chains.

Quick Start • Build & Usage

Installation, Compilation & Usage Examples

hydra_v2_sa is distributed as a self-contained standard C++17 source archive (pub/hydra_v2_sa.tgz) with zero external library dependencies beyond OpenMP and x86-64 AVX2/POPCNT intrinsics. Download the archive, extract it, enter the hydra_v2_sa directory, and run make to build the optimized ./hydra executable.

# Extract the archive, enter the directory, and compile the ./hydra binary

tar -xvzf hydra_v2_sa.tgz

cd hydra_v2_sa 
 
make



# Equivalent direct compiler command (inside hydra_v2_sa/):

g++ -O3 -mavx2 -mpopcnt -fopenmp -std=c++17 -Wall -Wextra src/hydra.cpp -o hydra

# 1. Default sequential Simulated Annealing on a benchmark alignment

./hydra examples/tst01.txt



# 2. Specify a custom PRNG seed and export the best topology in Newick format

./hydra -s 1 -o best_tree.nwk examples/ANGI.txt



# 3. Multi-chain Parallel SA on a 500-taxon alignment with 8 threads & 30s cap

./hydra --parallel -j 8 -t 30 -o zilla_best.nwk examples/zilla.txt



# 4. Custom annealing schedule in quiet scripting mode (prints integer score only)

./hydra -T 3.0 -c 0.95 -L 2500 -R 10 -q examples/CARP.txt

Supported Input Alignment Formats

PHYLIP-like (.txt): Optional header line <n_taxa> <n_chars> followed by whitespace-separated <taxon_name> <sequence> rows.
TNT / Hennig86 (.sg): Optional xread header and dimension line followed by taxon sequence rows (terminated optionally by ;).

Output & Score Verification

Exact Post-Check: Before reporting, HYDRA recomputes a full scratch bottom-up Fitch pass (compute_full_fitch()) to guarantee exact parsimony score integrity.
Newick Export (-o): Writes standard nested parenthetical tree topologies terminated by ; for downstream visualization (FigTree, iTOL, R ape).

Experimental Evaluation • Embedded results.txt snapshot

Benchmark Results: Sequential vs. Parallel Simulated Annealing

Tests were performed on an AMD Ryzen 5 9600X processor with 6 cores and 12 threads, 5.4 GHz, AVX2/AVX512

Comparative evaluation over 10 independent runs per instance (seeds 1..10) comparing default Non-Parallel (Sequential) SA against Parallel Multi-Chain SA (--parallel) across the 20 synthetic/biological benchmark instances (tst01.txt–tst20.txt, 45 to 75 taxa) and the large 500-taxon plant plastid rbcL benchmark (zilla.txt, 500 taxa × 759 sites).

Cumulative Min Score (tst01–20)

20,032 (−7)

--parallel vs. 20,039 sequential (improves min on 6/20)

Cumulative Avg Score (tst01–20)

20,036.70 (−38.20)

Improves mean parsimony on 19/20 instances; 10/20 hit 0 variance

Mean Runtime (tst01–20)

1.048 s

Sequential default (vs. 7.234 s for 4× oversubscribed --parallel)

Large Scale: zilla.txt (500 × 759)

16,218 100% hit

Sequential: 16,218.20 avg in 4.391 s • Parallel: 16,218.00 in 7.405 s

Benchmark Suite 1 — Standard Test Instances (`tst01.txt` to `tst20.txt`, 10 seeds each)

Show improved Min only

Instance (taxa × chars)	Non-Parallel / Sequential (Default)				Parallel OpenMP (`--parallel`)
Instance (taxa × chars)	Min	Average	Max	Avg Time	Min	Average	Max	Avg Time
tst01.txt 45×61	545	545.50	547	0.738 s	545	545.00 σ=0	545	5.037 s
tst02.txt 47×151	1354	1355.80	1358	0.878 s	1354	1354.10	1355	6.506 s
tst03.txt 49×111	833	834.00	835	0.863 s	833	833.00 σ=0	833	5.821 s
tst04.txt 50×97	588	588.10	589	0.886 s	587	587.80	588	6.263 s
tst05.txt 52×75	789	789.00	789	0.860 s	789	789.00 σ=0	789	5.795 s
tst06.txt 54×65	596	596.60	598	0.888 s	596	596.00 σ=0	596	5.926 s
tst07.txt 56×143	1269	1270.40	1272	1.078 s	1269	1269.00 σ=0	1269	7.626 s
tst08.txt 57×119	852	853.60	857	1.003 s	852	852.10	853	6.819 s
tst09.txt 59×93	1143	1144.20	1146	0.952 s	1141	1142.80	1143	6.702 s
tst10.txt 60×71	720	721.00	722	0.962 s	720	720.00 σ=0	720	6.668 s
tst11.txt 62×63	541	542.20	543	1.001 s	540	540.60	542	6.300 s
tst12.txt 64×147	1208	1214.60	1218	1.270 s	1208	1208.50	1211	8.689 s
tst13.txt 65×113	1515	1517.00	1519	1.061 s	1515	1515.00 σ=0	1515	7.733 s
tst14.txt 67×99	1161	1163.10	1167	1.100 s	1160	1160.20	1161	7.670 s
tst15.txt 69×77	753	754.40	756	1.150 s	752	752.00 σ=0	752	7.757 s
tst16.txt 70×69	527	530.30	533	1.108 s	527	527.20	529	7.473 s
tst17.txt 71×159	2450	2454.20	2460	1.384 s	2450	2450.00 σ=0	2450	10.087 s
tst18.txt 73×117	1521	1522.30	1526	1.292 s	1521	1521.00 σ=0	1521	8.802 s
tst19.txt 74×95	1013	1015.90	1018	1.240 s	1012	1012.30	1013	8.509 s
tst20.txt 75×79	661	662.70	664	1.247 s	661	661.10	662	8.486 s
Average Moyenne	1002.0	1003.75	1005.9	1.048 s	1001.6	1001.84	1002.4	7.234 s
Total 20 instances	20039	20074.90	20117	20.961 s	20032	20036.70	20047	144.671 s

Benchmark Suite 2 — Large-Scale Phylogenetic Dataset (`zilla.txt`, 500 Taxa × 759 Sites)

Greedy 1stRotuGbr Init (n ≥ 200)

Instance (taxa × chars)	Non-Parallel / Sequential (Default)				Parallel OpenMP (`--parallel`)
Instance (taxa × chars)	Min	Average	Max	Avg Time	Min	Average	Max	Avg Time
zilla.txt 500×759	16218	16218.20	16219	4.391 s	16218	16218.00 σ=0	16218	7.405 s

Analysis of Experimental Findings

Speed vs. Multi-Chain Precision Trade-off: Sequential SA solves the entire 20-instance tst01–tst20 suite in 20.961 s total (1.048 s mean per instance) while already reaching within 0.19% of the multi-chain parallel average.
Superior Solution Quality & Stability with --parallel: Launching 4 × hardware_threads cooperative Markov chains for n < 200 lowers cumulative minimum parsimony score from 20,039 to 20,032 (discovering strictly better best-known trees on tst04, tst09, tst11, tst14, tst15, and tst19), lowers cumulative average score by 38.20 steps (20,074.90 → 20,036.70), and achieves zero variance (Min = Avg = Max) on 11 out of 20 instances.
Scalability on zilla.txt (500 Taxa): Combining 1stRotuGbr greedy initialization with AVX2 two-pass SPR search reaches the optimal 16,218 parsimony score in 4.391 s sequentially (16,218.20 mean) and achieves a 100% hit rate (16,218.00 across all 10 seeds) in 7.405 s under --parallel (which scales at 1 × hardware_threads for n ≥ 200).

View Verbatim Benchmark File (results.txt) Click to expand / collapse

               |       Non parallèle (Séquentiel)        |         Parallèle (--parallel)         
Instance       |    Min |   Moyenne |    Max | Temps moy |    Min |   Moyenne |    Max | Temps moy
---------------+--------+-----------+--------+-----------+--------+-----------+--------+----------
tst01.txt      |    545 |    545.50 |    547 |   0.738 s |    545 |    545.00 |    545 |   5.037 s
tst02.txt      |   1354 |   1355.80 |   1358 |   0.878 s |   1354 |   1354.10 |   1355 |   6.506 s
tst03.txt      |    833 |    834.00 |    835 |   0.863 s |    833 |    833.00 |    833 |   5.821 s
tst04.txt      |    588 |    588.10 |    589 |   0.886 s |    587 |    587.80 |    588 |   6.263 s
tst05.txt      |    789 |    789.00 |    789 |   0.860 s |    789 |    789.00 |    789 |   5.795 s
tst06.txt      |    596 |    596.60 |    598 |   0.888 s |    596 |    596.00 |    596 |   5.926 s
tst07.txt      |   1269 |   1270.40 |   1272 |   1.078 s |   1269 |   1269.00 |   1269 |   7.626 s
tst08.txt      |    852 |    853.60 |    857 |   1.003 s |    852 |    852.10 |    853 |   6.819 s
tst09.txt      |   1143 |   1144.20 |   1146 |   0.952 s |   1141 |   1142.80 |   1143 |   6.702 s
tst10.txt      |    720 |    721.00 |    722 |   0.962 s |    720 |    720.00 |    720 |   6.668 s
tst11.txt      |    541 |    542.20 |    543 |   1.001 s |    540 |    540.60 |    542 |   6.300 s
tst12.txt      |   1208 |   1214.60 |   1218 |   1.270 s |   1208 |   1208.50 |   1211 |   8.689 s
tst13.txt      |   1515 |   1517.00 |   1519 |   1.061 s |   1515 |   1515.00 |   1515 |   7.733 s
tst14.txt      |   1161 |   1163.10 |   1167 |   1.100 s |   1160 |   1160.20 |   1161 |   7.670 s
tst15.txt      |    753 |    754.40 |    756 |   1.150 s |    752 |    752.00 |    752 |   7.757 s
tst16.txt      |    527 |    530.30 |    533 |   1.108 s |    527 |    527.20 |    529 |   7.473 s
tst17.txt      |   2450 |   2454.20 |   2460 |   1.384 s |   2450 |   2450.00 |   2450 |  10.087 s
tst18.txt      |   1521 |   1522.30 |   1526 |   1.292 s |   1521 |   1521.00 |   1521 |   8.802 s
tst19.txt      |   1013 |   1015.90 |   1018 |   1.240 s |   1012 |   1012.30 |   1013 |   8.509 s
tst20.txt      |    661 |    662.70 |    664 |   1.247 s |    661 |    661.10 |    662 |   8.486 s
---------------+--------+-----------+--------+-----------+--------+-----------+--------+----------
Moyenne        | 1002.0 |   1003.75 | 1005.9 |   1.048 s | 1001.6 |   1001.84 | 1002.4 |   7.234 s
Total          |  20039 |  20074.90 |  20117 |  20.961 s |  20032 |  20036.70 |  20047 | 144.671 s


               |       Non parallèle (Séquentiel)        |         Parallèle (--parallel)         
Instance       |    Min |   Moyenne |    Max | Temps moy |    Min |   Moyenne |    Max | Temps moy
---------------+--------+-----------+--------+-----------+--------+-----------+--------+----------
zilla.txt      |  16218 |  16218.20 |  16219 |   4.391 s |  16218 |  16218.00 |  16218 |   7.405 s

Bibliography

Scientific References

Richer, J.-M., Vazquez-Ortiz, K. E., & Lesaint, D. (2016). Strategies for Maximum Parsimony Phylogenetic Reconstruction. Electronic Notes in Discrete Mathematics (ENDM), 52, 165–172.
Richer, J.-M., Goëffon, A., & Hao, J.-K. (2008). Three New Techniques to Estimate a Maximum Parsimony Tree. (Introduces the randomized greedy 1stRotuGbr stepwise construction operator.)
Goëffon, A. (2006). Nouvelles heuristiques de voisinage et mémétique pour le problème de Maximum de Parcimonie. Ph.D. Thesis, Université d'Angers (Section 4.4.4).