nx_voronoi.nx source

1// nx_voronoi.nx -- 2D Voronoi diagram (procgen LAYER 3). 2// 3// Sites are scattered in [0, w) x [0, h); each cell of the canvas is 4// assigned to its NEAREST site by squared Euclidean distance. Output 5// is a per-pixel label array (which site index owns each pixel). 6// 7// USE CASES: 8// - City-block partition for procgen city generation (each block 9// becomes a Voronoi cell around its district seed) 10// - Zombie-escape agent territories (each spawn point owns the 11// region closest to it) 12// - Reaction-diffusion / cellular-automaton substrate (cell shape) 13// - Crystal-grain texture for procgen surfaces 14// - Lloyd relaxation seed for centroidal Voronoi tessellation 15// 16// ALGORITHM: 17// Brute-force nearest-site per pixel. For each pixel (x, y), scan 18// all n_sites and pick the index with minimum (dx^2 + dy^2). Cost: 19// O(w * h * n_sites). Patent-clean (squared-distance is textbook). 20// 21// SCALING NOTE: 22// Fortune's sweepline (1986) and Aurenhammer (1991) give O(n log n) 23// for the EDGES of the diagram, but per-pixel labelling for a 24// raster canvas is still O(w*h*n) without spatial acceleration. 25// For substrate use (canvas <= 256x256, sites <= 256), brute-force 26// is fast under qemu-riscv64 -- millions of i64 ops, no f64, no 27// sqrt. Upgrade path: jump-flooding-algorithm (JFA, Rong+Tan 2006) 28// for GPU-equivalent O(w*h*log(max(w,h))) when canvas grows. 29// 30// QUALITATIVE BAND (dual-reading cardinal): 31// Per-site cell-area sealed enum classifies how uniform the 32// partition is, derived from variance of per-site pixel count. 33// 34// COMPOSITION: 35// sites can be supplied by any source -- in practice 36// nx_poisson_disk_sample gives the most visually-pleasing 37// uniform-but-not-grid distribution. Demo wiring in 38// nx_voronoi_test.nx. 39// 40// genealogy_id: voronoi_1908 + dirichlet_1850 + fortune_1986 + 41// aurenhammer_1991 + rong_tan_2006_jfa 42// lineage_id: nx_voronoi_q10 43 44// nx_safety_envelope: 45// intended_use: AUTO_APPLIED -- primitive-specific tuning queued 46// sil_target: SIL1 47// evidence: [bulk_applied_2026-05-16, see-file-comment-for-detail] 48// verdict: NOT_YET_EVALUATED 49 50import "nx_syscalls.nx" 51import "nx_tier.nx" 52 53// ===== Sealed-enum: cell-area uniformity bands =================== 54// 55// After labelling, we compute per-site pixel counts and classify 56// the distribution's coefficient of variation (sd / mean) into one 57// of four bands. Q10 fixed-point throughout. 58 59const NX_VORONOI_BAND_UNIFORM: nx_int = 0 // very even, CV < 0.15 60const NX_VORONOI_BAND_BALANCED: nx_int = 1 // moderate spread, CV < 0.35 61const NX_VORONOI_BAND_LOPSIDED: nx_int = 2 // some sites dominate, CV < 0.65 62const NX_VORONOI_BAND_DEGENERATE: nx_int = 3 // wildly uneven, CV >= 0.65 63const NX_VORONOI_N_BANDS: nx_int = 4 64 65const NX_VORONOI_Q: nx_int = 1024 66// CV thresholds in Q10 (0.15, 0.35, 0.65 * 1024). 67const NX_VORONOI_CV_UNIFORM_Q10: nx_int = 154 68const NX_VORONOI_CV_BALANCED_Q10: nx_int = 358 69const NX_VORONOI_CV_LOPSIDED_Q10: nx_int = 666 70 71func nx_voronoi_band_is_valid(band: nx_int) -> nx_int { 72 if band < 0 { return 0 } 73 if band >= NX_VORONOI_N_BANDS { return 0 } 74 return 1 75} 76 77// ===== Per-pixel labelling ========================================= 78// 79// Writes labels[y * w + x] = nearest-site index for each pixel. 80// Pixels strictly outside (w, h) are not touched. Returns the 81// number of pixels labelled. 82 83func nx_voronoi_label(w: nx_int, h: nx_int, 84 site_xs: *i64, site_ys: *i64, n_sites: nx_int, 85 labels: *i64) -> nx_int { 86 if w <= 0 { return 0 } 87 if h <= 0 { return 0 } 88 if n_sites <= 0 { return 0 } 89 90 var y: nx_int = 0 91 while y < h { 92 var x: nx_int = 0 93 while x < w { 94 var best_idx: nx_int = 0 95 let dx0: nx_int = x - site_xs[0] 96 let dy0: nx_int = y - site_ys[0] 97 var best_d_sq: nx_int = dx0 * dx0 + dy0 * dy0 98 var i: nx_int = 1 99 while i < n_sites { 100 let dx: nx_int = x - site_xs[i] 101 let dy: nx_int = y - site_ys[i] 102 let d_sq: nx_int = dx * dx + dy * dy 103 if d_sq < best_d_sq { 104 best_d_sq = d_sq 105 best_idx = i 106 } 107 i = i + 1 108 } 109 labels[y * w + x] = best_idx 110 x = x + 1 111 } 112 y = y + 1 113 } 114 return w * h 115} 116 117// ===== Per-site pixel count ======================================= 118// 119// Walk the label grid once, accumulating per-site counts. Caller 120// supplies a zero-initialised counts[n_sites] buffer. 121 122func nx_voronoi_site_counts(w: nx_int, h: nx_int, 123 labels: *i64, n_sites: nx_int, 124 counts: *i64) -> nx_int { 125 if n_sites <= 0 { return 0 } 126 var i: nx_int = 0 127 while i < n_sites { 128 counts[i] = 0 129 i = i + 1 130 } 131 var y: nx_int = 0 132 while y < h { 133 var x: nx_int = 0 134 while x < w { 135 let lbl: nx_int = labels[y * w + x] 136 if lbl >= 0 { 137 if lbl < n_sites { 138 counts[lbl] = counts[lbl] + 1 139 } 140 } 141 x = x + 1 142 } 143 y = y + 1 144 } 145 return 0 146} 147 148// ===== Integer sqrt (Newton-Raphson, monotone) ===================== 149// 150// For CV computation: sqrt(variance). Newton-Raphson with initial 151// estimate via log2_floor; converges in ~6 iterations for inputs 152// up to 2^40. Returns floor(sqrt(n)) for non-negative n. 153 154func _vor_isqrt(n: nx_int) -> nx_int { 155 if n <= 0 { return 0 } 156 if n < 2 { return 1 } 157 // Initial estimate: half the bit-width. 158 var x: nx_int = n 159 var shift: nx_int = 0 160 while x > 0 { 161 x = x / 2 162 shift = shift + 1 163 } 164 var r: nx_int = 1 165 var i: nx_int = 0 166 while i < (shift / 2 + 1) { 167 r = r * 2 168 i = i + 1 169 } 170 // Newton iterations: r_next = (r + n/r) / 2 171 var iter: nx_int = 0 172 while iter < 20 { 173 if r <= 0 { return 0 } 174 let next: nx_int = (r + n / r) / 2 175 if next >= r { return r } 176 r = next 177 iter = iter + 1 178 } 179 return r 180} 181 182// ===== Coefficient-of-variation classifier ======================== 183// 184// Given counts[n_sites], compute mean, variance, sd, and CV (in Q10); 185// map CV into one of the four sealed bands. 186 187func nx_voronoi_classify(counts: *i64, n_sites: nx_int) -> nx_int { 188 if n_sites <= 0 { return NX_VORONOI_BAND_DEGENERATE } 189 var total: nx_int = 0 190 var i: nx_int = 0 191 while i < n_sites { 192 total = total + counts[i] 193 i = i + 1 194 } 195 if total <= 0 { return NX_VORONOI_BAND_DEGENERATE } 196 let mean: nx_int = total / n_sites 197 if mean <= 0 { return NX_VORONOI_BAND_DEGENERATE } 198 // Sum of squared deviations. 199 var sum_sq_dev: nx_int = 0 200 var j: nx_int = 0 201 while j < n_sites { 202 let dev: nx_int = counts[j] - mean 203 sum_sq_dev = sum_sq_dev + dev * dev 204 j = j + 1 205 } 206 let variance: nx_int = sum_sq_dev / n_sites 207 let sd: nx_int = _vor_isqrt(variance) 208 // CV in Q10 = (sd * Q) / mean. 209 let cv_q10: nx_int = (sd * NX_VORONOI_Q) / mean 210 if cv_q10 < NX_VORONOI_CV_UNIFORM_Q10 { return NX_VORONOI_BAND_UNIFORM } 211 if cv_q10 < NX_VORONOI_CV_BALANCED_Q10 { return NX_VORONOI_BAND_BALANCED } 212 if cv_q10 < NX_VORONOI_CV_LOPSIDED_Q10 { return NX_VORONOI_BAND_LOPSIDED } 213 return NX_VORONOI_BAND_DEGENERATE 214} 215 216// ===== Accessor: label lookup ===================================== 217// 218// Bounds-checked read; returns -1 if (x, y) is out of canvas. 219 220func nx_voronoi_label_at(labels: *i64, w: nx_int, h: nx_int, 221 x: nx_int, y: nx_int) -> nx_int { 222 if x < 0 { return 0 - 1 } 223 if x >= w { return 0 - 1 } 224 if y < 0 { return 0 - 1 } 225 if y >= h { return 0 - 1 } 226 return labels[y * w + x] 227}