708 lines
29 KiB
Rust
708 lines
29 KiB
Rust
//! Equilibrium hydrology solver benchmarks (T-1177, body-map-viewer workshop
|
||
//! measurement ①).
|
||
//!
|
||
//! Measures `hydrology_equilibrium::solve` at today's Layer-1 working grid
|
||
//! (512×256) and at two 4K-class synthetic grids (~768×432 ≈ 330K cells,
|
||
//! matching the workshop's per-gridunit derive measurement ②'s canvas size
|
||
//! for direct comparison; 3840×2160 ≈ 8.3M cells, the "computer catches
|
||
//! fire" ceiling case). Real GJ1c heightmap data is used at 512×256 (the
|
||
//! actual production working-grid size — no upsampling needed there); the
|
||
//! two larger grids use synthetic elevation (documented in
|
||
//! `synthetic_elevation` below) since no committed heightmap PNG is stored
|
||
//! at those resolutions and generating/committing new fixture PNGs is out of
|
||
//! scope for a measurement prototype.
|
||
//!
|
||
//! Run: `cargo test --release --test hydrology_equilibrium_bench -- --ignored --nocapture`
|
||
//! (debug numbers are not representative — this crate's other benches use
|
||
//! the same release-only convention).
|
||
//!
|
||
//! Hardware: 16 cores, Rayon default thread pool (14 workers observed
|
||
//! elsewhere in this repo's benches on the same machine).
|
||
|
||
use std::time::Instant;
|
||
|
||
use settled_reach_server::atlas::heightmap::load_heightmap_png;
|
||
use settled_reach_server::atlas::hydrology_equilibrium::{solve, ClimateInputs};
|
||
|
||
/// Deterministic synthetic elevation for grids larger than any committed
|
||
/// heightmap PNG. NOT a real body — a smooth multi-octave ridged surface
|
||
/// (a few sine terms at different frequencies/phases, summed and
|
||
/// normalized) chosen to produce a realistic MIX of basins or the solver
|
||
/// would have nothing to fill: a plain gradient (as `zoom_ladder_bench.rs`'s
|
||
/// `bench_hm` uses for its unrelated per-cell derive cost) has almost no
|
||
/// interior depressions, which would make this bench measure an
|
||
/// unrepresentative best case (priority-flood on a monotonic slope is
|
||
/// nearly free — the expensive part is basin interiors + overflow search).
|
||
/// Purely a function of `(row, col, width, height)` — the same call always
|
||
/// produces the same bytes, so the resulting elevation grid is itself
|
||
/// deterministic (D-010), even though it is synthetic rather than sourced
|
||
/// from a real body.
|
||
fn synthetic_elevation(width: u32, height: u32) -> Vec<f32> {
|
||
let w = width as f64;
|
||
let h = height as f64;
|
||
let n = (width * height) as usize;
|
||
(0..n)
|
||
.map(|i| {
|
||
let row = (i / width as usize) as f64;
|
||
let col = (i % width as usize) as f64;
|
||
let x = col / w;
|
||
let y = row / h;
|
||
// Several sine octaves at different frequencies/phases — enough
|
||
// basins (local minima not at the grid boundary) that the
|
||
// priority-flood + overflow-search work is representative, not
|
||
// a degenerate monotonic slope.
|
||
let v = 0.5
|
||
+ 0.25
|
||
* (x * std::f64::consts::TAU * 3.0).sin()
|
||
* (y * std::f64::consts::TAU * 2.0).cos()
|
||
+ 0.15
|
||
* (x * std::f64::consts::TAU * 7.3 + 1.7).sin()
|
||
* (y * std::f64::consts::TAU * 5.1).sin()
|
||
+ 0.10
|
||
* (x * std::f64::consts::TAU * 13.0).cos()
|
||
* (y * std::f64::consts::TAU * 11.0 + 0.4).sin();
|
||
v.clamp(0.0, 1.0) as f32
|
||
})
|
||
.collect()
|
||
}
|
||
|
||
fn gj1c_512x256() -> (Vec<f32>, f32) {
|
||
let src = std::path::PathBuf::from(env!("CARGO_MANIFEST_DIR"))
|
||
.join("../wiki/star-systems/GJ-1/bodies/GJ1c/heightmap.png");
|
||
let heightmap = load_heightmap_png(&src, "GJ1c", 0.3).expect("decode committed GJ1c heightmap");
|
||
let small = heightmap.downsample(512, 256); // GRID_W x GRID_H, the real production working grid
|
||
(small.data, small.sea_level)
|
||
}
|
||
|
||
fn default_climate() -> ClimateInputs {
|
||
ClimateInputs { moisture_q: 55 }
|
||
}
|
||
|
||
fn run_and_report(label: &str, width: u32, height: u32, elevation: &[f32], sea_level: f32) {
|
||
let n_cells = (width as u64) * (height as u64);
|
||
|
||
// Cold run.
|
||
let t0 = Instant::now();
|
||
let result_cold = solve(elevation, width, height, sea_level, default_climate());
|
||
let cold = t0.elapsed();
|
||
|
||
// Warm run (same process, allocator/cache warm — same input).
|
||
let t1 = Instant::now();
|
||
let result_warm = solve(elevation, width, height, sea_level, default_climate());
|
||
let warm = t1.elapsed();
|
||
|
||
let lake_cells: usize = result_cold.basins.iter().map(|b| b.cells.len()).sum();
|
||
let carved_cells = result_cold.cliff_edge.iter().filter(|&&c| c).count();
|
||
let endorheic_count = result_cold
|
||
.basins
|
||
.iter()
|
||
.filter(|b| {
|
||
!b.cells.is_empty()
|
||
&& matches!(
|
||
b.outcome,
|
||
settled_reach_server::atlas::hydrology_equilibrium::BasinOutcome::Endorheic { .. }
|
||
)
|
||
})
|
||
.count();
|
||
let overflow_count = result_cold
|
||
.basins
|
||
.iter()
|
||
.filter(|b| {
|
||
!b.cells.is_empty()
|
||
&& matches!(
|
||
b.outcome,
|
||
settled_reach_server::atlas::hydrology_equilibrium::BasinOutcome::Overflow { .. }
|
||
)
|
||
})
|
||
.count();
|
||
|
||
println!("\n=== {label} ({width}x{height} = {n_cells} cells) ===");
|
||
println!(
|
||
" cold: {:>9.2} ms total, {:>8.1} ns/cell",
|
||
cold.as_secs_f64() * 1000.0,
|
||
cold.as_secs_f64() * 1e9 / n_cells as f64
|
||
);
|
||
println!(
|
||
" warm: {:>9.2} ms total, {:>8.1} ns/cell",
|
||
warm.as_secs_f64() * 1000.0,
|
||
warm.as_secs_f64() * 1e9 / n_cells as f64
|
||
);
|
||
println!(
|
||
" basins: {} total ({} overflow, {} endorheic, {} empty/no-depression), \
|
||
lake cells: {lake_cells}, carved gorge cells: {carved_cells}",
|
||
result_cold.basins.len(),
|
||
overflow_count,
|
||
endorheic_count,
|
||
result_cold.basins.len() - overflow_count - endorheic_count,
|
||
);
|
||
std::hint::black_box(&result_warm);
|
||
}
|
||
|
||
#[test]
|
||
#[ignore]
|
||
fn bench_512x256_real_gj1c() {
|
||
let (elev, sea_level) = gj1c_512x256();
|
||
run_and_report(
|
||
"512x256 (real GJ1c, production working-grid size)",
|
||
512,
|
||
256,
|
||
&elev,
|
||
sea_level,
|
||
);
|
||
}
|
||
|
||
#[test]
|
||
#[ignore]
|
||
fn bench_768x432_synthetic() {
|
||
let (w, h) = (768u32, 432u32);
|
||
let elev = synthetic_elevation(w, h);
|
||
run_and_report(
|
||
"768x432 (~330K cells, 4K-class synthetic — see synthetic_elevation doc)",
|
||
w,
|
||
h,
|
||
&elev,
|
||
0.35,
|
||
);
|
||
}
|
||
|
||
#[test]
|
||
#[ignore]
|
||
fn bench_3840x2160_synthetic() {
|
||
let (w, h) = (3840u32, 2160u32);
|
||
let elev = synthetic_elevation(w, h);
|
||
run_and_report(
|
||
"3840x2160 (~8.3M cells, 4K synthetic — see synthetic_elevation doc)",
|
||
w,
|
||
h,
|
||
&elev,
|
||
0.35,
|
||
);
|
||
}
|
||
|
||
/// Determinism proof at bench scale (T-1177 mandatory deliverable): same
|
||
/// seed + input → byte-identical solver output, twice, on a non-trivial
|
||
/// grid (not just the small fixtures already covered by the module's own
|
||
/// unit tests).
|
||
#[test]
|
||
#[ignore]
|
||
fn determinism_at_330k_cells() {
|
||
let (w, h) = (768u32, 432u32);
|
||
let elev = synthetic_elevation(w, h);
|
||
let r1 = solve(&elev, w, h, 0.35, default_climate());
|
||
let r2 = solve(&elev, w, h, 0.35, default_climate());
|
||
assert_eq!(
|
||
r1.filled_scaled, r2.filled_scaled,
|
||
"filled surface must be byte-identical"
|
||
);
|
||
assert_eq!(
|
||
r1.channel_depth_scaled, r2.channel_depth_scaled,
|
||
"carved channel depth must be byte-identical"
|
||
);
|
||
assert_eq!(
|
||
r1.cliff_edge, r2.cliff_edge,
|
||
"cliff-edge flags must be byte-identical"
|
||
);
|
||
assert_eq!(
|
||
r1.basins.len(),
|
||
r2.basins.len(),
|
||
"basin count must be identical"
|
||
);
|
||
for (a, b) in r1.basins.iter().zip(r2.basins.iter()) {
|
||
assert_eq!(a.basin_id, b.basin_id);
|
||
assert_eq!(a.cells, b.cells);
|
||
assert_eq!(a.spill_level_scaled, b.spill_level_scaled);
|
||
assert_eq!(a.spill_cell, b.spill_cell);
|
||
assert_eq!(format!("{:?}", a.outcome), format!("{:?}", b.outcome));
|
||
}
|
||
println!(
|
||
"\n=== determinism proof (768x432, {} basins) — byte-identical across two solves ===",
|
||
r1.basins.len()
|
||
);
|
||
}
|
||
|
||
/// Rayon-parallel throughput: the REAL production shape is N independent
|
||
/// bodies, each solved once (not one body's solve parallelized internally —
|
||
/// priority-flood's heap and the overflow Dijkstra search are both globally
|
||
/// sequential by nature, same as `road_graph.rs`'s A*). This measures what
|
||
/// "always keep hydrology for ~273 bodies" would cost in wall-clock if
|
||
/// solved across the Rayon pool, at the 512×256 production grid size —
|
||
/// directly answering the workshop's red-flag-2-adjacent question of
|
||
/// whether per-body-open hydrology is affordable at scale.
|
||
///
|
||
/// **Scope note (population-survey gap, closed by
|
||
/// [`bench_population_survey_all_committed_bodies`] below):** this bench
|
||
/// solves ONE real body (GJ1c) 273 times — it proves the *wall-clock*
|
||
/// affordability of running 273 independent solves in parallel (the
|
||
/// question it was built to answer), but the 273 solves are not 273
|
||
/// *distinct* bodies, so it says nothing about how basin/endorheic/carved-
|
||
/// cell distributions vary across the real body population. That is a
|
||
/// different question, answered by the population survey, not this bench —
|
||
/// left as-is (not rewritten) since it still correctly answers the
|
||
/// question it was designed for.
|
||
#[test]
|
||
#[ignore]
|
||
fn bench_parallel_273_bodies_at_512x256() {
|
||
use rayon::prelude::*;
|
||
|
||
let (elev, sea_level) = gj1c_512x256();
|
||
let body_count = 273usize;
|
||
|
||
let t0 = Instant::now();
|
||
let total_basins: usize = (0..body_count)
|
||
.into_par_iter()
|
||
.map(|_| {
|
||
let result = solve(&elev, 512, 256, sea_level, default_climate());
|
||
result.basins.len()
|
||
})
|
||
.sum();
|
||
let elapsed = t0.elapsed();
|
||
|
||
println!(
|
||
"\n=== {body_count} bodies x 512x256, Rayon par_iter ({} threads available) ===",
|
||
std::thread::available_parallelism()
|
||
.map(|n| n.get())
|
||
.unwrap_or(0)
|
||
);
|
||
println!(
|
||
" {:>9.2} ms total, {:>7.2} ms/body average, {total_basins} basins summed",
|
||
elapsed.as_secs_f64() * 1000.0,
|
||
elapsed.as_secs_f64() * 1000.0 / body_count as f64
|
||
);
|
||
}
|
||
|
||
/// **POPULATION SURVEY (post-adversarial, Troblum finding).** Loads and
|
||
/// solves EVERY real committed heightmap PNG in the repo (267 distinct
|
||
/// bodies as of this bench, `wiki/star-systems/*/bodies/*/heightmap.png`),
|
||
/// each downsampled to the real 512×256 production working grid with its
|
||
/// own PNG-embedded `sea_level`, run independently in parallel across the
|
||
/// Rayon pool — the actual population-scale question
|
||
/// [`bench_parallel_273_bodies_at_512x256`] does not answer (that bench
|
||
/// solves ONE body 273 times, not 273 distinct bodies; see its doc comment
|
||
/// above). This closes the gap: total wall time for the REAL population,
|
||
/// and — the load-bearing part for the cliff Phase-4 ruling — the
|
||
/// distribution of basins/endorheic-basins/carved-cliff-edge-cells ACROSS
|
||
/// the population, not just on GJ1c.
|
||
///
|
||
/// **Known scope limit, stated explicitly (not hidden):** every body uses
|
||
/// the SAME `ClimateInputs { moisture_q: 55 }` (`default_climate()`) — the
|
||
/// module's own doc already states real per-body moisture requires wiring
|
||
/// in `BodyWorldState.districts`/`regions` climate context, which is out of
|
||
/// scope for this survey (same limitation the original T-1177 prototype
|
||
/// documented, "No per-basin moisture/climate lookup"). This means the
|
||
/// endorheic classification specifically is measured under a uniform
|
||
/// climate assumption, not each body's real hydrosphere/atmosphere-derived
|
||
/// moisture — a real per-body climate wiring pass could shift the
|
||
/// endorheic counts. It does NOT limit the carved-cliff-edge finding below:
|
||
/// carving depends on ELEVATION geometry (the narrow two-basin-saddle
|
||
/// condition documented in the module), not on `moisture_q` at all — the
|
||
/// climate input only gates the endorheic/overflow classification of an
|
||
/// already-identified basin, never whether carving occurs.
|
||
#[test]
|
||
#[ignore]
|
||
fn bench_population_survey_all_committed_bodies() {
|
||
use rayon::prelude::*;
|
||
use settled_reach_server::atlas::hydrology_equilibrium::{BasinOutcome, DownstreamTarget};
|
||
|
||
let wiki_root =
|
||
std::path::PathBuf::from(env!("CARGO_MANIFEST_DIR")).join("../wiki/star-systems");
|
||
|
||
// Discover every committed heightmap PNG (deterministic ordering: sort
|
||
// by path so the survey's own reporting order is stable run-to-run —
|
||
// not load-bearing for solve() itself, which is pure per-body, but
|
||
// keeps the printed output diffable).
|
||
let mut heightmap_paths: Vec<std::path::PathBuf> = Vec::new();
|
||
fn walk(dir: &std::path::Path, out: &mut Vec<std::path::PathBuf>) {
|
||
let Ok(entries) = std::fs::read_dir(dir) else {
|
||
return;
|
||
};
|
||
for entry in entries.flatten() {
|
||
let path = entry.path();
|
||
if path.is_dir() {
|
||
walk(&path, out);
|
||
} else if path.file_name().and_then(|n| n.to_str()) == Some("heightmap.png") {
|
||
out.push(path);
|
||
}
|
||
}
|
||
}
|
||
walk(&wiki_root, &mut heightmap_paths);
|
||
heightmap_paths.sort();
|
||
|
||
assert!(
|
||
heightmap_paths.len() > 200,
|
||
"expected the real committed body population (~267 as of this bench), found {} — \
|
||
did the wiki_root path resolve correctly? ({wiki_root:?})",
|
||
heightmap_paths.len()
|
||
);
|
||
|
||
struct BodySurvey {
|
||
body_id: String,
|
||
basins: usize,
|
||
overflow: usize,
|
||
endorheic: usize,
|
||
lake_cells: usize,
|
||
carved_cells: usize,
|
||
cliff_edge_cells: usize,
|
||
solve_ms: f64,
|
||
}
|
||
|
||
let t0 = Instant::now();
|
||
let surveys: Vec<BodySurvey> = heightmap_paths
|
||
.par_iter()
|
||
.map(|path| {
|
||
// body_id = the directory name one level up from heightmap.png
|
||
// (wiki/star-systems/<system>/bodies/<body_id>/heightmap.png).
|
||
let body_id = path
|
||
.parent()
|
||
.and_then(|p| p.file_name())
|
||
.and_then(|n| n.to_str())
|
||
.unwrap_or("UNKNOWN")
|
||
.to_string();
|
||
|
||
let heightmap =
|
||
settled_reach_server::atlas::heightmap::load_heightmap_png(path, &body_id, 0.3)
|
||
.unwrap_or_else(|e| {
|
||
panic!("decode committed heightmap for {body_id} ({path:?}): {e}")
|
||
});
|
||
let small = heightmap.downsample(512, 256);
|
||
|
||
let t_solve = Instant::now();
|
||
let result = solve(&small.data, 512, 256, small.sea_level, default_climate());
|
||
let solve_ms = t_solve.elapsed().as_secs_f64() * 1000.0;
|
||
|
||
let lake_cells: usize = result.basins.iter().map(|b| b.cells.len()).sum();
|
||
let carved_cells: usize = result
|
||
.basins
|
||
.iter()
|
||
.filter(|b| {
|
||
matches!(
|
||
b.outcome,
|
||
BasinOutcome::Overflow {
|
||
downstream_target: DownstreamTarget::Sea
|
||
| DownstreamTarget::Basin(_)
|
||
| DownstreamTarget::OpenSpillway,
|
||
..
|
||
}
|
||
) && result.channel_depth_scaled[b.spill_cell] > 0
|
||
})
|
||
.count();
|
||
let cliff_edge_cells = result.cliff_edge.iter().filter(|&&c| c).count();
|
||
let endorheic = result
|
||
.basins
|
||
.iter()
|
||
.filter(|b| matches!(b.outcome, BasinOutcome::Endorheic { .. }))
|
||
.count();
|
||
let overflow = result
|
||
.basins
|
||
.iter()
|
||
.filter(|b| matches!(b.outcome, BasinOutcome::Overflow { .. }))
|
||
.count();
|
||
|
||
BodySurvey {
|
||
body_id,
|
||
basins: result.basins.len(),
|
||
overflow,
|
||
endorheic,
|
||
lake_cells,
|
||
carved_cells,
|
||
cliff_edge_cells,
|
||
solve_ms,
|
||
}
|
||
})
|
||
.collect();
|
||
let elapsed = t0.elapsed();
|
||
|
||
let body_count = surveys.len();
|
||
let total_basins: usize = surveys.iter().map(|s| s.basins).sum();
|
||
let total_overflow: usize = surveys.iter().map(|s| s.overflow).sum();
|
||
let total_endorheic: usize = surveys.iter().map(|s| s.endorheic).sum();
|
||
let total_lake_cells: usize = surveys.iter().map(|s| s.lake_cells).sum();
|
||
let total_carved_cells: usize = surveys.iter().map(|s| s.carved_cells).sum();
|
||
let total_cliff_edge_cells: usize = surveys.iter().map(|s| s.cliff_edge_cells).sum();
|
||
let bodies_with_any_carving = surveys.iter().filter(|s| s.cliff_edge_cells > 0).count();
|
||
|
||
let mut by_cliff_edge: Vec<&BodySurvey> = surveys.iter().collect();
|
||
by_cliff_edge.sort_by(|a, b| b.cliff_edge_cells.cmp(&a.cliff_edge_cells));
|
||
|
||
println!(
|
||
"\n=== POPULATION SURVEY: {body_count} real committed bodies x 512x256, Rayon par_iter ({} threads available) ===",
|
||
std::thread::available_parallelism()
|
||
.map(|n| n.get())
|
||
.unwrap_or(0)
|
||
);
|
||
println!(
|
||
" {:>9.2} ms total wall time, {:>7.3} ms/body average",
|
||
elapsed.as_secs_f64() * 1000.0,
|
||
elapsed.as_secs_f64() * 1000.0 / body_count as f64
|
||
);
|
||
println!(
|
||
" basins: {total_basins} total across population ({total_overflow} overflow, {total_endorheic} endorheic)"
|
||
);
|
||
println!(" lake cells (summed across population): {total_lake_cells}");
|
||
println!(
|
||
" carved-outlet basins (channel_depth_scaled > 0 at spill cell, summed): {total_carved_cells}"
|
||
);
|
||
println!(" cliff_edge=true cells (summed across population): {total_cliff_edge_cells}");
|
||
println!(
|
||
" bodies with ANY cliff_edge cell: {bodies_with_any_carving} / {body_count} ({:.2}%)",
|
||
100.0 * bodies_with_any_carving as f64 / body_count as f64
|
||
);
|
||
println!("\n Top 15 bodies by cliff_edge cell count (the max-carve outliers):");
|
||
for s in by_cliff_edge.iter().take(15) {
|
||
println!(
|
||
" {:>10} cliff_edge_cells={:>4} carved_basins={:>2} basins={:>3} (overflow={:>3}, endorheic={:>2}) lake_cells={:>6} solve={:.2}ms",
|
||
s.body_id, s.cliff_edge_cells, s.carved_cells, s.basins, s.overflow, s.endorheic, s.lake_cells, s.solve_ms
|
||
);
|
||
}
|
||
let zero_carve_bodies = surveys.iter().filter(|s| s.cliff_edge_cells == 0).count();
|
||
println!(
|
||
"\n Bodies with ZERO cliff_edge cells: {zero_carve_bodies} / {body_count} ({:.2}%)",
|
||
100.0 * zero_carve_bodies as f64 / body_count as f64
|
||
);
|
||
|
||
// Determinism spot-check on the max-carve outlier (if any carving was
|
||
// found at all) — re-solve it once more and confirm byte-identical
|
||
// cliff_edge/channel_depth output, so the survey's headline outlier
|
||
// isn't itself an artifact of non-determinism.
|
||
if let Some(top) = by_cliff_edge.first() {
|
||
if top.cliff_edge_cells > 0 {
|
||
let path = heightmap_paths
|
||
.iter()
|
||
.find(|p| {
|
||
p.parent()
|
||
.and_then(|d| d.file_name())
|
||
.and_then(|n| n.to_str())
|
||
== Some(top.body_id.as_str())
|
||
})
|
||
.expect("outlier body path must exist (found via the same walk above)");
|
||
let heightmap =
|
||
settled_reach_server::atlas::heightmap::load_heightmap_png(path, &top.body_id, 0.3)
|
||
.expect("re-decode outlier heightmap for determinism spot-check");
|
||
let small = heightmap.downsample(512, 256);
|
||
let r1 = solve(&small.data, 512, 256, small.sea_level, default_climate());
|
||
let r2 = solve(&small.data, 512, 256, small.sea_level, default_climate());
|
||
assert_eq!(
|
||
r1.cliff_edge, r2.cliff_edge,
|
||
"max-carve outlier {} must be deterministic (cliff_edge)",
|
||
top.body_id
|
||
);
|
||
assert_eq!(
|
||
r1.channel_depth_scaled, r2.channel_depth_scaled,
|
||
"max-carve outlier {} must be deterministic (channel_depth_scaled)",
|
||
top.body_id
|
||
);
|
||
println!(
|
||
"\n Determinism spot-check on max-carve outlier ({}): PASS (byte-identical across two solves)",
|
||
top.body_id
|
||
);
|
||
}
|
||
}
|
||
}
|
||
|
||
/// **PER-BASIN SIZE DISTRIBUTION (body-map-viewer workshop, Araminta's
|
||
/// endorheic-carrier crossover question, routed via the coordinator).**
|
||
/// [`bench_population_survey_all_committed_bodies`] above reports only
|
||
/// AGGREGATE lake-cell totals (2,694,012 summed across 22,270 basins) — it
|
||
/// does not report individual basin SIZES, which is exactly what a
|
||
/// sparse-list-vs-dense-field crossover analysis needs (a sparse per-cell
|
||
/// list is cheap for small basins and expensive for large ones; the
|
||
/// crossover point depends on the actual size DISTRIBUTION, not the
|
||
/// aggregate). This bench solves every real committed body at the real
|
||
/// 512×256 production working grid (same population, same grid size as the
|
||
/// population survey above) and records every individual basin's cell
|
||
/// count, then reports the full distribution: histogram buckets, percentiles,
|
||
/// and the largest basins found, so the crossover call is made from a real
|
||
/// distribution rather than the single synthetic 8.3M-cell grid's one
|
||
/// data point (412,700 cells — T-1177's original measured appendix, a
|
||
/// SYNTHETIC ridged-terrain grid at a canvas size no production path
|
||
/// actually derives hydrology at synchronously, not a real-population
|
||
/// statistic).
|
||
#[test]
|
||
#[ignore]
|
||
fn bench_per_basin_size_distribution_real_population() {
|
||
use rayon::prelude::*;
|
||
|
||
let wiki_root =
|
||
std::path::PathBuf::from(env!("CARGO_MANIFEST_DIR")).join("../wiki/star-systems");
|
||
let mut heightmap_paths: Vec<std::path::PathBuf> = Vec::new();
|
||
fn walk(dir: &std::path::Path, out: &mut Vec<std::path::PathBuf>) {
|
||
let Ok(entries) = std::fs::read_dir(dir) else {
|
||
return;
|
||
};
|
||
for entry in entries.flatten() {
|
||
let path = entry.path();
|
||
if path.is_dir() {
|
||
walk(&path, out);
|
||
} else if path.file_name().and_then(|n| n.to_str()) == Some("heightmap.png") {
|
||
out.push(path);
|
||
}
|
||
}
|
||
}
|
||
walk(&wiki_root, &mut heightmap_paths);
|
||
heightmap_paths.sort();
|
||
assert!(
|
||
heightmap_paths.len() > 200,
|
||
"expected the real committed body population (~267), found {}",
|
||
heightmap_paths.len()
|
||
);
|
||
|
||
// (body_id, basin_cell_count, outcome_label) for EVERY basin across the
|
||
// whole population — the raw material for the distribution below.
|
||
let all_basins: Vec<(String, usize, &'static str)> = heightmap_paths
|
||
.par_iter()
|
||
.flat_map(|path| {
|
||
let body_id = path
|
||
.parent()
|
||
.and_then(|p| p.file_name())
|
||
.and_then(|n| n.to_str())
|
||
.unwrap_or("UNKNOWN")
|
||
.to_string();
|
||
let heightmap = settled_reach_server::atlas::heightmap::load_heightmap_png(
|
||
path, &body_id, 0.3,
|
||
)
|
||
.unwrap_or_else(|e| panic!("decode committed heightmap for {body_id} ({path:?}): {e}"));
|
||
let small = heightmap.downsample(512, 256);
|
||
let result = solve(&small.data, 512, 256, small.sea_level, default_climate());
|
||
result
|
||
.basins
|
||
.iter()
|
||
.map(|b| {
|
||
let outcome_label = match b.outcome {
|
||
settled_reach_server::atlas::hydrology_equilibrium::BasinOutcome::Overflow {
|
||
..
|
||
} => "overflow",
|
||
settled_reach_server::atlas::hydrology_equilibrium::BasinOutcome::Endorheic {
|
||
..
|
||
} => "endorheic",
|
||
};
|
||
(body_id.clone(), b.cells.len(), outcome_label)
|
||
})
|
||
.collect::<Vec<_>>()
|
||
})
|
||
.collect();
|
||
|
||
assert_eq!(
|
||
all_basins.len(),
|
||
22_270,
|
||
"basin count must match the population survey's own total exactly \
|
||
(same population, same grid size, same solver — a mismatch here \
|
||
means the two benches drifted, not that either is wrong)"
|
||
);
|
||
|
||
let mut sizes: Vec<usize> = all_basins.iter().map(|(_, n, _)| *n).collect();
|
||
sizes.sort_unstable();
|
||
|
||
let total: usize = sizes.iter().sum();
|
||
let n = sizes.len();
|
||
let percentile = |p: f64| -> usize {
|
||
let idx = ((p / 100.0) * (n - 1) as f64).round() as usize;
|
||
sizes[idx.min(n - 1)]
|
||
};
|
||
|
||
// Histogram buckets chosen to bracket plausible sparse-vs-dense
|
||
// crossover points: a Vec<(u16,u16)> cell entry is 4 raw bytes/cell
|
||
// (MessagePack-framed, somewhat less) vs. the dense `water` field's
|
||
// PNG-per-field rate of ~1.9 bytes/cell (T-1179) applied ONLY to that
|
||
// basin's own bounding extent — the buckets below span from "clearly
|
||
// sparse wins" (tiny basins) through "clearly dense wins" (the largest
|
||
// basins in this population).
|
||
let buckets: [(usize, usize, &str); 7] = [
|
||
(0, 10, "1-10 cells"),
|
||
(11, 50, "11-50 cells"),
|
||
(51, 200, "51-200 cells"),
|
||
(201, 1_000, "201-1,000 cells"),
|
||
(1_001, 5_000, "1,001-5,000 cells"),
|
||
(5_001, 20_000, "5,001-20,000 cells"),
|
||
(20_001, usize::MAX, "20,001+ cells"),
|
||
];
|
||
let mut bucket_counts = [0usize; 7];
|
||
let mut bucket_cell_totals = [0usize; 7];
|
||
for &size in &sizes {
|
||
for (i, (lo, hi, _)) in buckets.iter().enumerate() {
|
||
if size >= *lo && size <= *hi {
|
||
bucket_counts[i] += 1;
|
||
bucket_cell_totals[i] += size;
|
||
break;
|
||
}
|
||
}
|
||
}
|
||
|
||
let mut by_size: Vec<&(String, usize, &str)> = all_basins.iter().collect();
|
||
by_size.sort_by(|a, b| b.1.cmp(&a.1));
|
||
|
||
println!("\n=== PER-BASIN SIZE DISTRIBUTION: {n} basins, real 267-body population, 512x256 working grid ===");
|
||
println!(" total lake cells (cross-check vs population survey's 2,694,012): {total}");
|
||
println!(
|
||
" min={}, max={}, mean={:.1}, median(p50)={}",
|
||
sizes[0],
|
||
sizes[n - 1],
|
||
total as f64 / n as f64,
|
||
percentile(50.0)
|
||
);
|
||
println!(
|
||
" percentiles: p10={} p25={} p50={} p75={} p90={} p95={} p99={} p99.9={}",
|
||
percentile(10.0),
|
||
percentile(25.0),
|
||
percentile(50.0),
|
||
percentile(75.0),
|
||
percentile(90.0),
|
||
percentile(95.0),
|
||
percentile(99.0),
|
||
percentile(99.9)
|
||
);
|
||
|
||
println!("\n Histogram (basin cell-count buckets):");
|
||
for (i, (_, _, label)) in buckets.iter().enumerate() {
|
||
println!(
|
||
" {label:<18}: {:>6} basins ({:>5.2}% of basins), {:>9} cells total ({:>5.2}% of lake cells)",
|
||
bucket_counts[i],
|
||
100.0 * bucket_counts[i] as f64 / n as f64,
|
||
bucket_cell_totals[i],
|
||
100.0 * bucket_cell_totals[i] as f64 / total as f64
|
||
);
|
||
}
|
||
|
||
println!("\n Top 20 largest basins (real population, real bodies):");
|
||
for (body_id, size, outcome) in by_size.iter().take(20) {
|
||
println!(" {body_id:>12} cells={size:>7} outcome={outcome}");
|
||
}
|
||
|
||
// Bytes-if-sparse vs bytes-if-dense-per-basin-extent crossover, using
|
||
// T-1179's measured PNG-per-field rate (638,382 B / 331,776 cells) as
|
||
// the dense-encoding baseline, and 4 raw bytes/cell (u16,u16 pair,
|
||
// MessagePack framing not applied — a conservative/pessimistic sparse
|
||
// estimate, since MessagePack's compact array framing would cost less,
|
||
// per T-1179's own "6.00 not 7 bytes/cell" finding for a different
|
||
// field shape) for the sparse Vec<(u16,u16)> carrier.
|
||
const PNG_BYTES_PER_CELL: f64 = 638_382.0 / 331_776.0;
|
||
const SPARSE_BYTES_PER_CELL: f64 = 4.0;
|
||
let mut crossover_basin_count = 0usize;
|
||
let mut sparse_wins_cells = 0usize;
|
||
let mut dense_wins_cells = 0usize;
|
||
for &size in &sizes {
|
||
let sparse_cost = size as f64 * SPARSE_BYTES_PER_CELL;
|
||
let dense_cost = size as f64 * PNG_BYTES_PER_CELL;
|
||
if sparse_cost <= dense_cost {
|
||
crossover_basin_count += 1;
|
||
sparse_wins_cells += size;
|
||
} else {
|
||
dense_wins_cells += size;
|
||
}
|
||
}
|
||
println!(
|
||
"\n Sparse-vs-dense crossover (sparse={SPARSE_BYTES_PER_CELL} B/cell raw vs dense={:.3} B/cell PNG-per-field):",
|
||
PNG_BYTES_PER_CELL
|
||
);
|
||
println!(
|
||
" sparse cheaper for {crossover_basin_count} / {n} basins ({:.2}%), covering {sparse_wins_cells} cells ({:.2}% of lake cells)",
|
||
100.0 * crossover_basin_count as f64 / n as f64,
|
||
100.0 * sparse_wins_cells as f64 / total as f64
|
||
);
|
||
println!(
|
||
" dense cheaper for {} / {n} basins ({:.2}%), covering {dense_wins_cells} cells ({:.2}% of lake cells)",
|
||
n - crossover_basin_count,
|
||
100.0 * (n - crossover_basin_count) as f64 / n as f64,
|
||
100.0 * dense_wins_cells as f64 / total as f64
|
||
);
|
||
println!();
|
||
}
|