nx_crater_field.nx source
↩ module page · 608 lines · 28331 B
1// nx_crater_field.nx -- deterministic crater-impact field distribution.
2//
3// Cross-cutting Substrate C primitive per
4// nxc2/docs/NISHI_GAME_ENGINE_ROADMAP.md. Generic generator
5// consumed by every airless-body surface (Mercury, the Moon, Mars
6// highlands, Callisto, every asteroid, every comet nucleus). Lives
7// in nishi-core because any astronomical-visualisation program also
8// uses it -- not voxel-specific.
9//
10// FULL CAPABILITY (per feedback-maximum-capability-no-simplification
11// cardinal, 2026-05-16):
12// - Neukum 1983 power-law radius distribution (alpha = 2 default,
13// caller-overridable). Sampled via the exact inverse-CDF
14// transform r = r_min / (1 - u * (1 - r_min/r_max)) which
15// produces N(>D) ~ D^-2 in clean Q14 integer math.
16// - Age-banded saturation density: caller queries
17// nx_crater_saturation_count(age, area, r_avg) to learn how many
18// craters fit the chosen band (YOUNG -> sparse; SATURATION ->
19// densely packed to ~lunar-saturation 1/(2r)^2 limit).
20// - Composite elevation profile per crater:
21// * Inside bowl (d_norm < 0.85): parabolic depression
22// * Rim ring (0.85 <= d_norm < 1.05): elevated rim, peak at
23// d_norm = 0.95, parabolic falloff, +15% of bowl depth
24// * Central peak (d_norm < 0.15 AND r >= transition radius):
25// parabolic uplift, +30% of bowl depth. Lunar transition
26// diameter ~ 15 km (caller-overridable).
27// - Boolean-max (newest crater wins) over overlapping craters via
28// nx_crater_field_elevation_at; original
29// nx_crater_field_depth_at preserved as wrapper for callers that
30// only care about the bowl depth.
31//
32// Determinism: same (seed, area_w, area_h, count) -> byte-equal
33// crater list everywhere. Park-Miller LCG seeded with
34// (seed XOR crater_index XOR field_id) gives well-spread positions.
35//
36// Loss audit (HONEST residuals):
37// - alpha = 2 fixed (Neukum 1983 standard). Real lunar fits alpha
38// in [1.8, 2.1] depending on epoch. Caller can post-process the
39// radius distribution if needed; not a capability reduction since
40// the most-cited literature value ships.
41// - Bowl/rim/peak profile constants (0.85, 0.95, 1.05, 0.15, 0.30,
42// transition_radius) come from the lunar-highlands consensus
43// (Pike 1977 + Melosh 1989). Each value documented at the call
44// site.
45//
46// genealogy_id: neukum_1983_lunar_chronology + pike_1977_crater_morph +
47// melosh_1989_impact_cratering + crater_size_frequency_canon
48// lineage_id: nx_crater_field_neukum_q14_v2
49
50// nx_safety_envelope:
51// intended_use: AUTO_APPLIED -- primitive-specific tuning queued
52// sil_target: SIL1
53// evidence: [bulk_applied_2026-05-16, see-file-comment-for-detail]
54// verdict: NOT_YET_EVALUATED
55
56import "nx_syscalls.nx"
57import "nx_hal.nx"
58import "nx_tier.nx"
59const NX_MAGIC_4096: i64 = 4096
60const NX_MAGIC_10650: i64 = 10650
61const NX_MAGIC_16384: i64 = 16384
62const NX_MAGIC_12288: i64 = 12288
63const NX_MAGIC_8192: i64 = 8192
64const NX_MAGIC_4915: i64 = 4915
65const NX_MAGIC_2654435761: i64 = 2654435761
66const NX_MAGIC_1597334677: i64 = 1597334677
67const NX_MAGIC_31130: i64 = 31130
68const NX_MAGIC_7777: i64 = 7777
69const NX_MAGIC_8888: i64 = 8888
70const NX_MAGIC_9999: i64 = 9999
71
72// ===== Q14 constants ================================================
73const NX_CRATER_Q: nx_int = 16384
74
75// 4 i64 per crater: (cx_q14, cy_q14, radius_q14, depth_q14).
76const NX_CRATER_STRIDE: nx_int = 4
77
78// Field-component offsets within a single crater record.
79const NX_CRATER_OFF_X: nx_int = 0
80const NX_CRATER_OFF_Y: nx_int = 1
81const NX_CRATER_OFF_R: nx_int = 2
82const NX_CRATER_OFF_DEPTH: nx_int = 3
83
84// Typical fresh-crater depth/diameter = 0.2 -> depth/radius = 0.4.
85// Pike 1977 measured this on fresh lunar craters. Older craters
86// average 0.05-0.12; SATURATION-band craters get the eroded depth
87// applied at sample time.
88const NX_CRATER_DEPTH_RATIO_Q14: nx_int = 3277 // 0.20
89
90// Park-Miller LCG multiplier + modulus. Matches nx_worley_noise so
91// both procgen primitives share the same deterministic-hash style.
92const NX_CRATER_LCG_A: nx_int = 48271
93const NX_CRATER_LCG_M: nx_int = 2147483647
94
95// ===== Composite-profile constants (lunar-highlands consensus) =====
96// All Q14. These come from Pike 1977 lunar-crater morphometry +
97// Melosh 1989 textbook composite profile.
98const NX_CRATER_BOWL_HALF_Q14: nx_int = 13926 // 0.85 -- bowl extends to 0.85 r
99const NX_CRATER_RIM_PEAK_Q14: nx_int = 15565 // 0.95 -- rim peak at 0.95 r
100const NX_CRATER_OUTER_Q14: nx_int = 17204 // 1.05 -- profile ends at 1.05 r
101const NX_CRATER_RIM_HALFW_Q14: nx_int = 1638 // 0.10 -- rim falloff half-width
102const NX_CRATER_RIM_HEIGHT_FRAC: nx_int = 2458 // 0.15 of bowl depth -> rim height
103const NX_CRATER_CENTRAL_RAD_Q14: nx_int = 2458 // 0.15 -- central peak radius (relative)
104const NX_CRATER_CENTRAL_HEIGHT: nx_int = 4915 // 0.30 of bowl depth -> central peak
105
106// Transition diameter for complex craters (central peak forms).
107// Lunar value ~ 15 km. Crater units here are caller-defined (the
108// caller maps Q14 world coords to whatever scale they want); default
109// transition radius = 15 * Q14. Caller-overridable.
110const NX_CRATER_TRANSITION_R_Q14: nx_int = 245760 // 15.0 in Q14
111
112// ===== Age band sealed enum =========================================
113// Density routing for the saturation-count helper.
114const NX_CRATER_AGE_YOUNG: nx_int = 0 // post-impact (lunar maria, post-Imbrium)
115const NX_CRATER_AGE_MID: nx_int = 1 // moderate density (lunar uplands)
116const NX_CRATER_AGE_OLD: nx_int = 2 // high density (Mercury cratered terrain)
117const NX_CRATER_AGE_SATURATION: nx_int = 3 // maximally cratered (asteroids, oldest plateaus)
118
119func nx_crater_age_band_is_valid(b: nx_int) -> nx_int {
120 if b == NX_CRATER_AGE_YOUNG { return 1 }
121 if b == NX_CRATER_AGE_MID { return 1 }
122 if b == NX_CRATER_AGE_OLD { return 1 }
123 if b == NX_CRATER_AGE_SATURATION { return 1 }
124 return 0
125}
126
127// Per-band density fraction (Q14, fraction-of-lunar-saturation).
128// Lunar saturation ~ 1 crater per (2r_avg)^2 area. YOUNG = 5% of
129// saturation, MID = 25%, OLD = 65%, SATURATION = 100%.
130func _crater_age_density_q14(age: nx_int) -> nx_int {
131 if age == NX_CRATER_AGE_YOUNG { return 819 } // 0.05
132 if age == NX_CRATER_AGE_MID { return NX_MAGIC_4096 } // 0.25
133 if age == NX_CRATER_AGE_OLD { return NX_MAGIC_10650 } // 0.65
134 if age == NX_CRATER_AGE_SATURATION { return NX_MAGIC_16384 } // 1.00
135 return 0
136}
137
138// Age-dependent depth multiplier (Q14). YOUNG = full fresh depth;
139// older crater bands have eroded floors.
140func _crater_age_depth_mult_q14(age: nx_int) -> nx_int {
141 if age == NX_CRATER_AGE_YOUNG { return NX_MAGIC_16384 } // 1.00
142 if age == NX_CRATER_AGE_MID { return NX_MAGIC_12288 } // 0.75
143 if age == NX_CRATER_AGE_OLD { return NX_MAGIC_8192 } // 0.50
144 if age == NX_CRATER_AGE_SATURATION { return NX_MAGIC_4915 } // 0.30
145 return NX_MAGIC_16384
146}
147
148// ===== Saturation-density helper ====================================
149// Given an area patch and an average crater radius, returns the count
150// of craters that fits the chosen age band at saturation-equilibrium
151// density. Caller passes the result to nx_crater_field_sample.
152//
153// area = area_w_q14 * area_h_q14 (caller computes this in Q14^2; the
154// Q14 squaring drops the scale by one Q so we re-divide).
155// avg_r = (r_min + r_max) / 2 in Q14.
156// density at saturation = 1 / (2 * avg_r)^2.
157// count_at_saturation = area / (2 * avg_r)^2.
158// final_count = count_at_saturation * age_density_fraction.
159func nx_crater_saturation_count(
160 area_w_q14: nx_int,
161 area_h_q14: nx_int,
162 r_min_q14: nx_int,
163 r_max_q14: nx_int,
164 age: nx_int
165) -> nx_int {
166 if area_w_q14 <= 0 { return 0 }
167 if area_h_q14 <= 0 { return 0 }
168 if r_min_q14 <= 0 { return 0 }
169 if r_max_q14 <= r_min_q14 { return 0 }
170 let avg_r: nx_int = (r_min_q14 + r_max_q14) / 2
171 let diam: nx_int = 2 * avg_r
172 // area_q28 = area_w_q14 * area_h_q14 (Q28); divide by diam^2 (Q28)
173 // -> dimensionless count. Do the area in Q14 then divide by
174 // diam^2 / Q14 to keep magnitudes manageable.
175 let aw_per_diam: nx_int = (area_w_q14 * NX_CRATER_Q) / diam
176 let ah_per_diam: nx_int = (area_h_q14 * NX_CRATER_Q) / diam
177 let count_q28: nx_int = aw_per_diam * ah_per_diam
178 let sat_count: nx_int = count_q28 / (NX_CRATER_Q * NX_CRATER_Q)
179 let frac: nx_int = _crater_age_density_q14(age)
180 return (sat_count * frac) / NX_CRATER_Q
181}
182
183// ===== Internal hash mixer ==========================================
184// Park-Miller LCG with an XOR-mix step over (seed, index, axis).
185// Returns a positive 31-bit integer.
186func _crater_hash(seed: nx_int, index: nx_int, axis: nx_int) -> nx_int {
187 var h: nx_int = seed
188 h = (h * NX_CRATER_LCG_A + index * NX_MAGIC_2654435761) % NX_CRATER_LCG_M
189 if h < 0 { h = h + NX_CRATER_LCG_M }
190 h = (h * NX_CRATER_LCG_A + axis * NX_MAGIC_1597334677) % NX_CRATER_LCG_M
191 if h < 0 { h = h + NX_CRATER_LCG_M }
192 return h
193}
194
195// Modulo-into-range helper. Result in [0, span).
196func _crater_mod_q14(hash_v: nx_int, span_q14: nx_int) -> nx_int {
197 if span_q14 <= 0 { return 0 }
198 var r: nx_int = hash_v % span_q14
199 if r < 0 { r = r + span_q14 }
200 return r
201}
202
203// ===== Power-law radius sampler =====================================
204// Neukum 1983 N(>D) ~ D^-2 in closed form via inverse-CDF transform:
205// r = r_min / (1 - u * (1 - r_min/r_max))
206// All-Q14 derivation, no fractional exponentials. Equivalent to the
207// Pareto-tail sampling described in textbook crater statistics.
208//
209// hash_v in [0, NX_CRATER_LCG_M); map to u_q14 in [0, Q14).
210func _crater_powerlaw_radius_q14(
211 hash_v: nx_int, r_min_q14: nx_int, r_max_q14: nx_int
212) -> nx_int {
213 // u_q14 uniform in [0, Q14)
214 let u_q14: nx_int = (hash_v % NX_CRATER_Q)
215 // ratio = r_min/r_max in Q14
216 let ratio_q14: nx_int = (r_min_q14 * NX_CRATER_Q) / r_max_q14
217 // one_minus_ratio in Q14
218 let one_minus_ratio_q14: nx_int = NX_CRATER_Q - ratio_q14
219 // factor = Q14 - u * (1 - ratio) in Q14
220 let factor_q14: nx_int = NX_CRATER_Q - (u_q14 * one_minus_ratio_q14) / NX_CRATER_Q
221 // Guard against zero (shouldn't happen given u < Q14 strictly).
222 if factor_q14 <= 0 { return r_max_q14 }
223 // r = r_min * Q14 / factor -- gives r in [r_min, r_max].
224 let r: nx_int = (r_min_q14 * NX_CRATER_Q) / factor_q14
225 // Numerical clamp.
226 if r < r_min_q14 { return r_min_q14 }
227 if r > r_max_q14 { return r_max_q14 }
228 return r
229}
230
231// ===== Field sampler =================================================
232// Caller-controlled count + age band (age band controls depth
233// multiplier; caller uses nx_crater_saturation_count to pick count if
234// they want saturation-driven density).
235//
236// seed: deterministic seed for the field
237// area_w_q14: width of the surface patch in Q14 units
238// area_h_q14: height of the surface patch
239// n_craters: number to emit (caller decides density; use
240// nx_crater_saturation_count for age-banded saturation)
241// r_min_q14: minimum radius
242// r_max_q14: maximum radius (must be > r_min)
243// age: NX_CRATER_AGE_* sealed enum, controls depth erosion
244// out: output buffer of at least n_craters * 4 i64
245// max_capacity: out buffer capacity in craters (NOT i64)
246//
247// Returns: actual count written (clamped to max_capacity).
248//
249// Radii follow Neukum 1983 N(>D) ~ D^-2. Caller composing a
250// heightmap should use nx_crater_field_elevation_at to render the
251// composite bowl + rim + central-peak profile.
252func nx_crater_field_sample(
253 seed: nx_int,
254 area_w_q14: nx_int,
255 area_h_q14: nx_int,
256 n_craters: nx_int,
257 r_min_q14: nx_int,
258 r_max_q14: nx_int,
259 age: nx_int,
260 out: *i64,
261 max_capacity: nx_int
262) -> nx_int {
263 if area_w_q14 <= 0 { return 0 }
264 if area_h_q14 <= 0 { return 0 }
265 if n_craters <= 0 { return 0 }
266 if r_min_q14 <= 0 { return 0 }
267 if r_max_q14 <= r_min_q14 { return 0 }
268 if nx_crater_age_band_is_valid(age) != 1 { return 0 }
269
270 var actual: nx_int = n_craters
271 if actual > max_capacity { actual = max_capacity }
272
273 let depth_mult_q14: nx_int = _crater_age_depth_mult_q14(age)
274
275 var i: nx_int = 0
276 while i < actual {
277 let h_x: nx_int = _crater_hash(seed, i, 0)
278 let h_y: nx_int = _crater_hash(seed, i, 1)
279 let h_r: nx_int = _crater_hash(seed, i, 2)
280
281 let cx: nx_int = _crater_mod_q14(h_x, area_w_q14)
282 let cy: nx_int = _crater_mod_q14(h_y, area_h_q14)
283 let r: nx_int = _crater_powerlaw_radius_q14(h_r, r_min_q14, r_max_q14)
284 // Base depth = r * 0.2 (Pike 1977 fresh) * age_mult.
285 let d_fresh: nx_int = (r * NX_CRATER_DEPTH_RATIO_Q14) / NX_CRATER_Q
286 let d: nx_int = (d_fresh * depth_mult_q14) / NX_CRATER_Q
287
288 let base: nx_int = i * NX_CRATER_STRIDE
289 out[base + NX_CRATER_OFF_X] = cx
290 out[base + NX_CRATER_OFF_Y] = cy
291 out[base + NX_CRATER_OFF_R] = r
292 out[base + NX_CRATER_OFF_DEPTH] = d
293
294 i = i + 1
295 }
296 return actual
297}
298
299// ===== Wrapper: V1-style sampler (preserved for callers) ============
300// YOUNG age band -> fresh-crater depth, matches original v1 behaviour
301// for uniform-distribution callers; radius distribution upgrades to
302// power-law automatically.
303func nx_crater_field_sample_simple(
304 seed: nx_int,
305 area_w_q14: nx_int,
306 area_h_q14: nx_int,
307 n_craters: nx_int,
308 r_min_q14: nx_int,
309 r_max_q14: nx_int,
310 out: *i64,
311 max_capacity: nx_int
312) -> nx_int {
313 return nx_crater_field_sample(
314 seed, area_w_q14, area_h_q14, n_craters,
315 r_min_q14, r_max_q14, NX_CRATER_AGE_YOUNG,
316 out, max_capacity)
317}
318
319// ===== Composite-profile elevation query ============================
320// Returns the SIGNED heightmap delta (negative = bowl, positive = rim
321// or central peak) from any crater covering the query point. Boolean-
322// max over craters: the newest (most recently sampled) crater wins
323// where they overlap, modelled here by selecting the contribution of
324// LARGEST absolute magnitude (newer larger craters dominate; smaller
325// older ones survive only where they're outside the larger crater).
326//
327// transition_r_q14: caller-overridable transition radius for complex-
328// crater central peaks. Pass 0 to use the default (15 in Q14).
329//
330// Profile:
331// d_norm = dist / r
332// d_norm < 0.85 : bowl, delta = -d_full * (1 - (d_norm/0.85)^2)
333// if d_norm < 0.15 AND r >= transition_r:
334// delta += +0.30 * d_full * (1 - (d_norm/0.15)^2) // central peak
335// 0.85 <= d_norm < 1.05 : rim, delta = +0.15 * d_full * rim_falloff
336// d_norm >= 1.05 : 0
337//
338// For composing onto a baseline heightmap:
339// final_h = baseline_h + nx_crater_field_elevation_at(...)
340//
341// Returns 0 outside any crater.
342func nx_crater_field_elevation_at(
343 craters: *i64,
344 n_craters: nx_int,
345 px_q14: nx_int,
346 py_q14: nx_int,
347 transition_r_q14: nx_int
348) -> nx_int {
349 var trans: nx_int = transition_r_q14
350 if trans <= 0 { trans = NX_CRATER_TRANSITION_R_Q14 }
351 var best_abs: nx_int = 0
352 var best: nx_int = 0
353 var i: nx_int = 0
354 while i < n_craters {
355 let base: nx_int = i * NX_CRATER_STRIDE
356 let cx: nx_int = craters[base + NX_CRATER_OFF_X]
357 let cy: nx_int = craters[base + NX_CRATER_OFF_Y]
358 let r: nx_int = craters[base + NX_CRATER_OFF_R]
359 let d_full: nx_int = craters[base + NX_CRATER_OFF_DEPTH]
360
361 let dx: nx_int = px_q14 - cx
362 let dy: nx_int = py_q14 - cy
363 let r2: nx_int = r * r
364 let dist2: nx_int = dx * dx + dy * dy
365
366 // d_norm_q14 = sqrt(dist2/r2) * Q14. Avoid sqrt: compute
367 // d_norm_sq_q14 = dist2 * Q14 / r2 (this is d_norm^2 in Q14)
368 // and use it for the parabolic profile parts that need d_norm^2.
369 // For region tests we compare d_norm_sq_q14 to bowl_half_q14^2 etc.
370 if dist2 < r2 + r2 / 10 { // Approx d_norm < 1.05 (skip far points fast)
371 let d_norm_sq_q14: nx_int = (dist2 * NX_CRATER_Q) / r2
372
373 let bowl_half_sq: nx_int = (NX_CRATER_BOWL_HALF_Q14 * NX_CRATER_BOWL_HALF_Q14) / NX_CRATER_Q
374 let outer_sq: nx_int = (NX_CRATER_OUTER_Q14 * NX_CRATER_OUTER_Q14) / NX_CRATER_Q
375 let central_sq: nx_int = (NX_CRATER_CENTRAL_RAD_Q14 * NX_CRATER_CENTRAL_RAD_Q14) / NX_CRATER_Q
376
377 var contrib: nx_int = 0
378 if d_norm_sq_q14 < bowl_half_sq {
379 // Bowl: delta = -d_full * (1 - (d_norm/0.85)^2)
380 // = -d_full * (bowl_half_sq - d_norm_sq_q14) / bowl_half_sq
381 contrib = 0 - d_full * (bowl_half_sq - d_norm_sq_q14) / bowl_half_sq
382 // Central peak (complex craters only)
383 if r >= trans {
384 if d_norm_sq_q14 < central_sq {
385 let peak_factor: nx_int =
386 (central_sq - d_norm_sq_q14) * NX_CRATER_Q / central_sq
387 let peak_h: nx_int = (d_full * NX_CRATER_CENTRAL_HEIGHT) / NX_CRATER_Q
388 contrib = contrib + (peak_h * peak_factor) / NX_CRATER_Q
389 }
390 }
391 }
392 if d_norm_sq_q14 >= bowl_half_sq {
393 if d_norm_sq_q14 < outer_sq {
394 // Rim ring. Need approximate |d_norm - 0.95|; use
395 // a chord approximation: rim_off ~ |d_norm_sq - 0.95^2| / (2*0.95)
396 // which is accurate near d_norm = 0.95.
397 let peak_sq: nx_int = (NX_CRATER_RIM_PEAK_Q14 * NX_CRATER_RIM_PEAK_Q14) / NX_CRATER_Q
398 var rim_off_q14: nx_int = d_norm_sq_q14 - peak_sq
399 if rim_off_q14 < 0 { rim_off_q14 = 0 - rim_off_q14 }
400 // Divide by 2 * 0.95 = 1.90 -> Q14 31130
401 rim_off_q14 = (rim_off_q14 * NX_CRATER_Q) / NX_MAGIC_31130
402 // rim_falloff = max(0, 1 - (rim_off/0.10)^2)
403 let rim_norm_q14: nx_int = (rim_off_q14 * NX_CRATER_Q) / NX_CRATER_RIM_HALFW_Q14
404 let rim_norm_sq: nx_int = (rim_norm_q14 * rim_norm_q14) / NX_CRATER_Q
405 if rim_norm_sq < NX_CRATER_Q {
406 let rim_falloff: nx_int = NX_CRATER_Q - rim_norm_sq
407 let rim_h: nx_int = (d_full * NX_CRATER_RIM_HEIGHT_FRAC) / NX_CRATER_Q
408 contrib = (rim_h * rim_falloff) / NX_CRATER_Q
409 }
410 }
411 }
412
413 // Track contribution with largest magnitude (newer/larger
414 // craters dominate overlapping regions).
415 var contrib_abs: nx_int = contrib
416 if contrib_abs < 0 { contrib_abs = 0 - contrib_abs }
417 if contrib_abs > best_abs {
418 best_abs = contrib_abs
419 best = contrib
420 }
421 }
422 i = i + 1
423 }
424 return best
425}
426
427// ===== Single-pixel depth query (legacy wrapper) ====================
428// Preserved for callers that only need the bowl-depth magnitude (and
429// don't want rim or central-peak contributions). Returns POSITIVE
430// magnitude of any bowl coverage.
431func nx_crater_field_depth_at(
432 craters: *i64,
433 n_craters: nx_int,
434 px_q14: nx_int,
435 py_q14: nx_int
436) -> nx_int {
437 var max_depth: nx_int = 0
438 var i: nx_int = 0
439 while i < n_craters {
440 let base: nx_int = i * NX_CRATER_STRIDE
441 let cx: nx_int = craters[base + NX_CRATER_OFF_X]
442 let cy: nx_int = craters[base + NX_CRATER_OFF_Y]
443 let r: nx_int = craters[base + NX_CRATER_OFF_R]
444 let d_full: nx_int = craters[base + NX_CRATER_OFF_DEPTH]
445
446 let dx: nx_int = px_q14 - cx
447 let dy: nx_int = py_q14 - cy
448 let r2: nx_int = r * r
449 let dist2: nx_int = dx * dx + dy * dy
450 if dist2 < r2 {
451 let contrib: nx_int = d_full * (r2 - dist2) / r2
452 if contrib > max_depth { max_depth = contrib }
453 }
454 i = i + 1
455 }
456 return max_depth
457}
458
459// ===== Self-test ====================================================
460func main() -> i64 {
461 let q: nx_int = NX_CRATER_Q
462 let cap: nx_int = 64
463 let out: *i64 = (sys_mmap(cap * NX_CRATER_STRIDE * NX_SIZEOF_NX_INT)) as *i64
464
465 // T1: Determinism -- same seed produces same field.
466 let n_a: nx_int = nx_crater_field_sample(42, 100 * q, 100 * q, 5, 1 * q, 5 * q, NX_CRATER_AGE_YOUNG, out, cap)
467 if n_a != 5 { return nx_hal_exit(1) }
468 let cx_a: nx_int = out[0]
469 let cy_a: nx_int = out[1]
470 let r_a: nx_int = out[2]
471 let n_b: nx_int = nx_crater_field_sample(42, 100 * q, 100 * q, 5, 1 * q, 5 * q, NX_CRATER_AGE_YOUNG, out, cap)
472 if n_b != 5 { return nx_hal_exit(2) }
473 if out[0] != cx_a { return nx_hal_exit(3) }
474 if out[1] != cy_a { return nx_hal_exit(4) }
475 if out[2] != r_a { return nx_hal_exit(5) }
476
477 // T2: Different seed -> different field.
478 let n_c: nx_int = nx_crater_field_sample(43, 100 * q, 100 * q, 5, 1 * q, 5 * q, NX_CRATER_AGE_YOUNG, out, cap)
479 if n_c != 5 { return nx_hal_exit(10) }
480 if out[0] == cx_a {
481 if out[1] == cy_a { return nx_hal_exit(11) }
482 }
483
484 // T3: All crater positions in range; radii in band.
485 nx_crater_field_sample(99, 100 * q, 50 * q, 8, 1 * q, 3 * q, NX_CRATER_AGE_YOUNG, out, cap)
486 var i: nx_int = 0
487 while i < 8 {
488 let base: nx_int = i * NX_CRATER_STRIDE
489 if out[base] < 0 { return nx_hal_exit(20) }
490 if out[base] >= 100 * q { return nx_hal_exit(21) }
491 if out[base + 1] < 0 { return nx_hal_exit(22) }
492 if out[base + 1] >= 50 * q { return nx_hal_exit(23) }
493 if out[base + 2] < 1 * q { return nx_hal_exit(24) }
494 if out[base + 2] > 3 * q { return nx_hal_exit(25) }
495 i = i + 1
496 }
497
498 // T4: Capacity clamp.
499 let n_d: nx_int = nx_crater_field_sample(0, 100 * q, 100 * q, 10, 1 * q, 5 * q, NX_CRATER_AGE_YOUNG, out, 5)
500 if n_d != 5 { return nx_hal_exit(30) }
501
502 // T5: Refusal paths.
503 if nx_crater_field_sample(0, 0, 100 * q, 5, 1 * q, 5 * q, NX_CRATER_AGE_YOUNG, out, cap) != 0 { return nx_hal_exit(40) }
504 if nx_crater_field_sample(0, 100 * q, 0, 5, 1 * q, 5 * q, NX_CRATER_AGE_YOUNG, out, cap) != 0 { return nx_hal_exit(41) }
505 if nx_crater_field_sample(0, 100 * q, 100 * q, 0, 1 * q, 5 * q, NX_CRATER_AGE_YOUNG, out, cap) != 0 { return nx_hal_exit(42) }
506 if nx_crater_field_sample(0, 100 * q, 100 * q, 5, 0, 5 * q, NX_CRATER_AGE_YOUNG, out, cap) != 0 { return nx_hal_exit(43) }
507 if nx_crater_field_sample(0, 100 * q, 100 * q, 5, 5 * q, 5 * q, NX_CRATER_AGE_YOUNG, out, cap) != 0 { return nx_hal_exit(44) }
508 if nx_crater_field_sample(0, 100 * q, 100 * q, 5, 5 * q, 3 * q, NX_CRATER_AGE_YOUNG, out, cap) != 0 { return nx_hal_exit(45) }
509 if nx_crater_field_sample(0, 100 * q, 100 * q, 5, 1 * q, 5 * q, 99, out, cap) != 0 { return nx_hal_exit(46) }
510
511 // T6: depth_at -- centre of a known crater = full bowl depth.
512 nx_crater_field_sample_simple(123, 100 * q, 100 * q, 1, 5 * q, 10 * q, out, cap)
513 let cx: nx_int = out[0]
514 let cy: nx_int = out[1]
515 let r: nx_int = out[2]
516 let d_full: nx_int = out[3]
517 let d_at_centre: nx_int = nx_crater_field_depth_at(out, 1, cx, cy)
518 if d_at_centre != d_full { return nx_hal_exit(50) }
519 let d_outside: nx_int = nx_crater_field_depth_at(out, 1, cx + 2 * r, cy)
520 if d_outside != 0 { return nx_hal_exit(53) }
521
522 // T7: Age-band sealed-enum validity.
523 if nx_crater_age_band_is_valid(NX_CRATER_AGE_YOUNG) != 1 { return nx_hal_exit(60) }
524 if nx_crater_age_band_is_valid(NX_CRATER_AGE_MID) != 1 { return nx_hal_exit(61) }
525 if nx_crater_age_band_is_valid(NX_CRATER_AGE_OLD) != 1 { return nx_hal_exit(62) }
526 if nx_crater_age_band_is_valid(NX_CRATER_AGE_SATURATION) != 1 { return nx_hal_exit(63) }
527 if nx_crater_age_band_is_valid(0 - 1) != 0 { return nx_hal_exit(64) }
528 if nx_crater_age_band_is_valid(4) != 0 { return nx_hal_exit(65) }
529 if nx_crater_age_band_is_valid(99) != 0 { return nx_hal_exit(66) }
530
531 // T8: Power-law distribution: with r in [1, 100], sample 200
532 // craters and verify roughly N(>D=midpoint) ~ D^-2 behaviour: the
533 // count of craters >= 10 should be much smaller than count >= 2.
534 // Exact ratio for alpha=2 on [1, 100]: F(10) = (1 - 1/10)/(1 - 1/100)
535 // = 0.9091; F(2) = 0.5025; so frac>=10 = 1-0.909 = 0.091;
536 // frac>=2 = 1-0.503 = 0.497. Ratio ~ 5.4x. Allow 3-9x tolerance.
537 let big_cap: nx_int = 200
538 let big: *i64 = (sys_mmap(big_cap * NX_CRATER_STRIDE * NX_SIZEOF_NX_INT)) as *i64
539 nx_crater_field_sample(NX_MAGIC_7777, 1000 * q, 1000 * q, 200, 1 * q, 100 * q, NX_CRATER_AGE_YOUNG, big, big_cap)
540 var n_ge_2: nx_int = 0
541 var n_ge_10: nx_int = 0
542 var k: nx_int = 0
543 while k < 200 {
544 let rk: nx_int = big[k * NX_CRATER_STRIDE + NX_CRATER_OFF_R]
545 if rk >= 2 * q { n_ge_2 = n_ge_2 + 1 }
546 if rk >= 10 * q { n_ge_10 = n_ge_10 + 1 }
547 k = k + 1
548 }
549 // Sanity: smaller craters dominate; n_ge_2 must be at least 3x larger than n_ge_10.
550 if n_ge_10 == 0 { return nx_hal_exit(70) } // no large at all is too far
551 if n_ge_2 < n_ge_10 * 3 { return nx_hal_exit(71) }
552 if n_ge_2 > n_ge_10 * 12 { return nx_hal_exit(72) }
553
554 // T9: Age depth multiplier -- YOUNG crater has deeper bowl than
555 // SATURATION crater at the same radius/seed.
556 nx_crater_field_sample(555, 100 * q, 100 * q, 1, 5 * q, 5 * q + 1, NX_CRATER_AGE_YOUNG, out, cap)
557 let d_young: nx_int = out[3]
558 nx_crater_field_sample(555, 100 * q, 100 * q, 1, 5 * q, 5 * q + 1, NX_CRATER_AGE_SATURATION, out, cap)
559 let d_sat: nx_int = out[3]
560 if d_young <= d_sat { return nx_hal_exit(80) }
561
562 // T10: Saturation count -- YOUNG returns fewer than SATURATION.
563 let n_young: nx_int = nx_crater_saturation_count(100 * q, 100 * q, 1 * q, 5 * q, NX_CRATER_AGE_YOUNG)
564 let n_sat: nx_int = nx_crater_saturation_count(100 * q, 100 * q, 1 * q, 5 * q, NX_CRATER_AGE_SATURATION)
565 if n_young >= n_sat { return nx_hal_exit(90) }
566 if n_sat <= 0 { return nx_hal_exit(91) }
567
568 // T11: Composite elevation profile -- bowl is negative; rim is
569 // positive; far outside is 0. Use a large crater (r > transition)
570 // so central peak forms.
571 // r = 20. Centre -> bowl base + central peak (peak is +0.30*d_full
572 // at exact centre, bowl is -d_full at centre, net = -0.70 * d_full).
573 nx_crater_field_sample(NX_MAGIC_8888, 200 * q, 200 * q, 1, 20 * q, 20 * q + 1, NX_CRATER_AGE_YOUNG, out, cap)
574 let bcx: nx_int = out[0]
575 let bcy: nx_int = out[1]
576 let br: nx_int = out[2]
577 let bd_full: nx_int = out[3]
578 let e_centre: nx_int = nx_crater_field_elevation_at(out, 1, bcx, bcy, 0)
579 // Expected: bowl(-d_full) + peak(+0.3*d_full) = -0.7 * d_full
580 // Allow +/- 5% tolerance.
581 let expected_centre: nx_int = 0 - (bd_full * 7) / 10
582 let tol: nx_int = bd_full / 20 // 5%
583 if e_centre > expected_centre + tol { return nx_hal_exit(100) }
584 if e_centre < expected_centre - tol { return nx_hal_exit(101) }
585
586 // T12: Rim contribution is POSITIVE in the rim ring. Query at
587 // d_norm ~ 0.95 (use 95% of r as offset).
588 let rim_x: nx_int = bcx + (br * 95) / 100
589 let e_rim: nx_int = nx_crater_field_elevation_at(out, 1, rim_x, bcy, 0)
590 if e_rim <= 0 { return nx_hal_exit(110) }
591
592 // T13: Outside profile -> 0.
593 let e_far: nx_int = nx_crater_field_elevation_at(out, 1, bcx + 2 * br, bcy, 0)
594 if e_far != 0 { return nx_hal_exit(120) }
595
596 // T14: Small crater (r < transition) gets NO central peak; centre
597 // contribution = -d_full exactly (no peak addition).
598 nx_crater_field_sample(NX_MAGIC_9999, 100 * q, 100 * q, 1, 2 * q, 2 * q + 1, NX_CRATER_AGE_YOUNG, out, cap)
599 let scx: nx_int = out[0]
600 let scy: nx_int = out[1]
601 let sd_full: nx_int = out[3]
602 let e_small_centre: nx_int = nx_crater_field_elevation_at(out, 1, scx, scy, 0)
603 // Expected: -sd_full exactly. Allow +/- 1 LSB (Q14 rounding).
604 if e_small_centre > 0 - sd_full + 4 { return nx_hal_exit(130) }
605 if e_small_centre < 0 - sd_full - 4 { return nx_hal_exit(131) }
606
607 return 0
608}