diff --git a/build.md b/build.md index 54c1e93..d6cfffb 100644 --- a/build.md +++ b/build.md @@ -35,11 +35,11 @@ make -j With no explicit `SPACETIME` setting, this builds all supported spacetimes: -| Executable | Spacetime | -| --- | --- | -| `build/Release/minkowski_sky` | Flat Minkowski spacetime | +| Executable | Spacetime | +| --------------------------------- | --------------------------------------------------------- | +| `build/Release/minkowski_sky` | Flat Minkowski spacetime | | `build/Release/schwarzschild_sky` | Analytic Schwarzschild in ingoing Kerr–Schild coordinates | -| `build/Release/alcubierre_sky` | Analytic moving Alcubierre warp bubble, `x_s(t)=v_s t` (no capture) | +| `build/Release/alcubierre_sky` | Analytic moving Alcubierre warp bubble, `x_s(t)=v_s t` | To build only one: diff --git a/nr_spacetime_movie_renderer_design.md b/nr_spacetime_movie_renderer_design.md index 616471d..a2e85f2 100644 --- a/nr_spacetime_movie_renderer_design.md +++ b/nr_spacetime_movie_renderer_design.md @@ -1130,8 +1130,6 @@ residual、Chebyshev 表或解析主项；运行期不得建表。 各 trial 先确定可表示的目标时间，再用实际 `target-before.t` 推进状态、构造 dense interpolant 与提交时间；不得以舍入前的步长推进空间、舍入后的步长记时间。 worldtube 采样的坐标时间同样使用实际差值回推运动，不假设半步总是可表示。 - Alcubierre 的默认时间范围由 `1.25*4*R_escape/(1-|v_s|)` 给出；固定步的步数 - 估算不能截短自适应路径的历史，也不能替代其独立资源上限。 - 初始 accepted 状态的 RHS/metric 失败直接报告具体 point reason，不缩步；后续 stage 的 `OUT_OF_DOMAIN/INVALID_METRIC` 可有界缩步重试；`TIME_UNAVAILABLE/INTERNAL_ERROR` 不靠无限缩步；最小步长或时间不可进导致失败时报告 `INCOMPLETE/INTEGRATION_ERROR`。 diff --git a/src/main.c b/src/main.c index 9fa029f..33291df 100644 --- a/src/main.c +++ b/src/main.c @@ -846,8 +846,8 @@ static void print_help(const char *program) { stdout); #ifdef SPACETIME_ALCUBIERRE fputs( - "\nAlcubierre warp bubble (moving x_s(t)=v_s*t; no capture):\n" - " --alcubierre-vs V Constant bubble velocity v_s, |v_s| < 1 (default: 0.5)\n" + "\nAlcubierre warp bubble (moving x_s(t)=v_s*t; no geometric capture, shared dark policy):\n" + " --alcubierre-vs V Constant bubble velocity v_s (any finite value, default: 0.5)\n" " --alcubierre-radius R Bubble radius R > 0 (default: 5)\n" " --alcubierre-sigma S Wall sharpness sigma > 0 (default: 1)\n" " The escape radius R + 20/sigma is derived internally; the\n" @@ -991,19 +991,57 @@ static double alcubierre_time_step(const Settings *s) { return fmin(0.1, 0.05 / s->alcubierre_sigma); } -/* Coordinate-time coverage a past-directed Alcubierre ray must be granted. - * A ray that is nearly comoving with the bubble separates from its center in - * the propagation direction at only ~1 - |v_s|, so crossing the ~4*escape - * domain takes up to ~4*escape/(1 - |v_s|), scaled by ALCUBIERRE_BUDGET_MARGIN - * for lingering in the wall. This is a physical budget derived from the bubble - * geometry and the asymptotic separation rate; it is independent of the chosen - * step size and of the accepted-step count. */ +/* Coordinate-time resource allowance for a past-directed Alcubierre ray: + * + * B = margin * 4 * R_escape / max(|1 - |v_s||, exp(-D)), D = dark threshold. + * + * For a sub-luminal ray whose separation |1 - |v_s|| dominates exp(-D) this + * reduces to the historical margin*4*escape/(1-|v_s|) geometry budget + * (margin*4 == 5). For + * |v_s| -> 1, including exactly |v_s| = 1 and super-luminal |v_s| > 1, the + * exp(-D) floor is motivated by the center-comoving axial relation + * exp(-(L-L0)) = 1 - |v_s|*(1-f); it cuts off the separation scale at the + * finite dark threshold. This is a resource heuristic, not a universal bound + * for arbitrary cameras or directions. B is a resource allowance, not a + * guarantee that every ray escapes: a ray that needs more coordinate time ends + * as UNRESOLVED/BUDGET_EXHAUSTED, which the render-level publication gate + * handles. B is independent of the chosen step size and accepted-step count. + * + * Saturation rules: + * - exp(-D) may underflow to 0 for a large D; the max still selects the + * finite separation unless the separation itself is 0. + * - when the ordinary 5*R/denom is not finite and positive (denom == 0, or + * the division/scale overflows), fall back to a log-space evaluation + * log B = log(5) + log(R) - max(log(separation), -D), which never forms + * exp(D); clamp the log to [log(DBL_MIN), log(DBL_MAX/4)] and exponentiate. + * - the ordinary result is clamped to the same [DBL_MIN, DBL_MAX/4] range, + * reserving DBL_MAX/4 for the default 4x retry growth. + * - exp(D) is never evaluated. */ static double alcubierre_trace_time_budget(const Settings *s) { const double escape = spacetime_alcubierre_escape_radius(s->alcubierre_radius, s->alcubierre_sigma); - const double separation = 1.0 - fabs(s->alcubierre_vs); - return ALCUBIERRE_BUDGET_MARGIN * 4.0 * escape / separation; + const double separation = fabs(1.0 - fabs(s->alcubierre_vs)); + const double floor = exp(-s->dark_threshold); /* may underflow to 0 */ + const double upper = DBL_MAX / 4.0; + if (!isfinite(escape) || escape <= 0.0) + return DBL_MIN; + const double denom = fmax(separation, floor); + if (denom > 0.0 && isfinite(denom)) { + /* Do not pre-scale 5*R (which could overflow); divide first. */ + const double scaled = (escape / denom) * (ALCUBIERRE_BUDGET_MARGIN * 4.0); + if (isfinite(scaled) && scaled > 0.0) + return fmin(fmax(scaled, DBL_MIN), upper); + } + /* Extreme fallback: separation == 0 with exp(-D) underflowed to 0, or the + * ordinary arithmetic saturated. Evaluate in log space without exp(D). */ + double log_budget = log(ALCUBIERRE_BUDGET_MARGIN * 4.0) + log(escape) - + fmax(log(separation), -s->dark_threshold); + if (!isfinite(log_budget)) + log_budget = log(upper); + log_budget = fmin(log_budget, log(upper)); + log_budget = fmax(log_budget, log(DBL_MIN)); + return fmin(fmax(exp(log_budget), DBL_MIN), upper); } /* Legacy accepted-step estimate for the same coverage, used only by the @@ -1190,10 +1228,13 @@ static GeodesicTraceConfig trace_config(const Settings *s) { #elif defined(SPACETIME_ALCUBIERRE) default_step = alcubierre_time_step(s); { - /* Accepted-step resource estimate along the legacy fixed-step policy. It - * may be capped at ALCUBIERRE_MAX_TRACE_STEPS; the coordinate-time coverage - * below is what guarantees the physical budget, so a capped count never - * shortens the trusted history. */ + /* Accepted-step resource estimate along the legacy fixed-step policy. + * `alcubierre_step_budget` may be +inf or overflow for a near-luminal or + * super-luminal bubble, so the estimate is only cast when it is finite and + * strictly below the cap; otherwise the cap is kept. A capped count is an + * independent resource allowance and never shortens the trusted history, + * which the coordinate-time coverage below defines. It does not guarantee + * escape for every ray. */ const double estimated_steps = alcubierre_step_budget(s); unsigned int budget = ALCUBIERRE_MAX_TRACE_STEPS; if (isfinite(estimated_steps) && estimated_steps < (double)budget) @@ -2393,33 +2434,28 @@ int main(int argc, char **argv) { if (spacetime_create_alcubierre(&spacetime, settings.alcubierre_vs, settings.alcubierre_radius, settings.alcubierre_sigma)) { - fputs("Could not create Alcubierre spacetime source; require |v_s| < 1, " - "R > 0, sigma > 0.\n", stderr); + fputs("Could not create Alcubierre spacetime source; require finite " + "v_s, R > 0, sigma > 0, and a derived escape radius beyond R.\n", + stderr); return 1; } /* The fixed-step RK4 comparison path sizes its history from - * step * accepted-step count, so its domain guard still rejects an - * estimate that exceeds the cap. The adaptive path uses an independent - * coordinate-time budget and must not be rejected here: trustworthy quota - * exhaustion is reported as UNRESOLVED and handled by the publication - * gate, and explicit lookback/step overrides remain independent. */ + * step * accepted-step count, so its derived default estimate is still + * bounded by the cap. An explicit --trace-max-steps is a user resource + * allowance and is never rejected here; the adaptive path uses an + * independent coordinate-time budget and is never rejected here either + * (trustworthy quota exhaustion is reported as UNRESOLVED and handled by + * the publication gate). */ if (settings.stepper == GEODESIC_STEPPER_RK4 && + !settings.trace_max_steps_specified && alcubierre_step_budget(&settings) > (double)ALCUBIERRE_MAX_TRACE_STEPS) { - /* The estimate has a V-shaped minimum at sigma = 0.5, where the step - * stops being capped: below it the 20/sigma term dominates (increase - * sigma helps), above it the step scales as 1/sigma (decrease sigma - * helps), and at exactly 0.5 neither direction improves anything. */ - const char *sigma_advice = ""; - if (settings.alcubierre_sigma > 0.5) - sigma_advice = "decrease --alcubierre-sigma, "; - else if (settings.alcubierre_sigma < 0.5) - sigma_advice = "increase --alcubierre-sigma, "; fprintf(stderr, - "Fixed-step RK4 Alcubierre trace budget exceeds the %u-step " - "cap; decrease --alcubierre-radius, %sor move --alcubierre-vs " - "away from +/-1, or use the adaptive default.\n", - ALCUBIERRE_MAX_TRACE_STEPS, sigma_advice); + "Fixed-step RK4 Alcubierre default trace budget exceeds the " + "%u-step cap; decrease --alcubierre-radius, adjust " + "--dark-threshold or --alcubierre-sigma, or use the adaptive " + "default or an explicit --trace-max-steps.\n", + ALCUBIERRE_MAX_TRACE_STEPS); spacetime_destroy(&spacetime); return 2; } diff --git a/src/spacetime.h b/src/spacetime.h index 34ad702..cf919f4 100644 --- a/src/spacetime.h +++ b/src/spacetime.h @@ -122,8 +122,10 @@ int spacetime_create_default(SpacetimeSource *source); int spacetime_create_minkowski(SpacetimeSource *source, double escape_radius); int spacetime_create_schwarzschild_ks(SpacetimeSource *source, double mass, double escape_radius); -/* Moving Alcubierre bubble with x_s(t) = vs*t and x_s(0) = 0. Requires - * |vs| < 1, R > 0, and sigma > 0. */ +/* Moving Alcubierre bubble with x_s(t) = vs*t and x_s(0) = 0. Requires a + * finite vs, R > 0, and sigma > 0. Sub- and super-luminal |vs| are accepted; + * classify() only ever reports ACTIVE or ESCAPED, and the shared + * camera-relative dark policy may terminate a ray as DARK for any finite vs. */ int spacetime_create_alcubierre(SpacetimeSource *source, double vs, double radius, double sigma); /* Bubble-centered escape radius used by the Alcubierre backend; also lets diff --git a/src/spacetime_alcubierre.c b/src/spacetime_alcubierre.c index 6a75be2..0376eba 100644 --- a/src/spacetime_alcubierre.c +++ b/src/spacetime_alcubierre.c @@ -107,11 +107,8 @@ static SpacetimePointStatus alcubierre_eval(const SpacetimeSource *source, return SPACETIME_POINT_OK; } -/* A warp bubble has no curvature singularity or horizon for |v_s| < 1, so - * rays are only ever ACTIVE or ESCAPED; the exotic matter that would source - * the bubble is treated as optically transparent. The escape sphere follows - * the bubble, so rays terminate only once the metric is flat to machine - * precision at their current location. */ +/* The exotic matter that would source the bubble is treated as optically + * transparent. The escape sphere follows the moving bubble. */ static SpacetimeRayStatus alcubierre_classify(const SpacetimeSource *source, double t, const double x[3]) { const AlcubierreContext *context = source->context; @@ -179,7 +176,7 @@ double spacetime_alcubierre_escape_radius(double radius, double sigma) { int spacetime_create_alcubierre(SpacetimeSource *source, double vs, double radius, double sigma) { - if (source == NULL || !isfinite(vs) || fabs(vs) >= 1.0 || + if (source == NULL || !isfinite(vs) || !isfinite(radius) || radius <= 0.0 || !isfinite(sigma) || sigma <= 0.0) return -1; /* Reject parameter combinations whose derived domain overflows or does not diff --git a/tests/test_adaptive_cli.py b/tests/test_adaptive_cli.py index c40718d..b632c69 100644 --- a/tests/test_adaptive_cli.py +++ b/tests/test_adaptive_cli.py @@ -11,6 +11,7 @@ fixed-step reference convergence check. Small CPU 16x8/32x16 scenes keep the runtime short; physical comparisons use stored lens-map endpoints, not only the rendered PNG. """ +import math import os import re import struct @@ -139,6 +140,28 @@ def sum_rhs(vertices): return sum(v[14] for v in vertices) +def alcubierre_budget(escape, vs, dark_threshold): + """Reference implementation of the production Alcubierre time allowance: + B = 5*escape / max(|1-|v_s||, exp(-D)), with the DBL_MIN..DBL_MAX/4 + saturation and the log-space fallback used when the ordinary division is + not finite and positive. exp(D) is never formed.""" + sep = abs(1.0 - abs(vs)) + floor = math.exp(-dark_threshold) + denom = max(sep, floor) + upper = sys.float_info.max / 4.0 + if denom > 0.0 and math.isfinite(denom): + scaled = 5.0 * (escape / denom) + if math.isfinite(scaled) and scaled > 0.0: + return min(max(scaled, sys.float_info.min), upper) + log_sep = math.log(sep) if sep > 0.0 else -math.inf + log_budget = math.log(5.0) + math.log(escape) - max(log_sep, -dark_threshold) + if not math.isfinite(log_budget): + log_budget = math.log(upper) + log_budget = min(log_budget, math.log(upper)) + log_budget = max(log_budget, math.log(sys.float_info.min)) + return min(max(math.exp(log_budget), sys.float_info.min), upper) + + with tempfile.TemporaryDirectory(prefix='gr-adaptive-cli-', dir=str(TMP_ROOT)) as directory: tmp = Path(directory) @@ -435,111 +458,217 @@ with tempfile.TemporaryDirectory(prefix='gr-adaptive-cli-', print(f'{backend}: adaptive CLI checks passed', flush=True) - # Alcubierre production-default smoke: the DP54 default policy must - # validate and the migrated lookback budget must actually cover the warp - # bubble feature (rays integrate and escape) instead of pre-routing all - # misses. No movie/thread sweep for this third backend. + # Alcubierre production policy: sub- and super-luminal velocities share one + # finite resource allowance B = 5*escape / max(|1-|v_s||, exp(-D)); the DP54 + # default must validate, actually integrate the warp feature, and reach both + # the shared DARK terminal and escapes. All images are tiny 8x4/16x8. alc = BUILD / 'alcubierre_sky' if not alc.exists(): print('alcubierre: binary absent, skipping', flush=True) else: alc_help = run(alc, '--help').stdout ext = 'png' if '.png' in alc_help else 'ppm' - alc_map = tmp / 'alcubierre_default.grlens' - alc_out = tmp / f'alcubierre_default.{ext}' - run(alc, '--alcubierre-vs', 0.3, '--alcubierre-radius', 1, - '--alcubierre-sigma', 1, '--catalog', 'assets/sky_grid_5deg.csv', - '--width', 16, '--height', 8, '--fov-deg', 80, '--exposure', 1e-3, - '--coarse-cell-pixels', 8, '--refine-max-level', 0, - '--psf-relative-tail', 1e-4, '--verbose', '--output', alc_out, - '--lens-map-output', alc_map) - _, _, alc_prov = map_provenance(alc_map) - assert alc_prov['integrator'] == 1, alc_prov - assert alc_prov['min_step'] <= alc_prov['coordinate_time_step'] \ - <= alc_prov['max_step'] - assert alc_prov['max_lookback_time'] > 0 - alc_vertices, _ = map_vertices(alc_map) - outcomes = {v[10] for v in alc_vertices} - assert 3 not in outcomes, f'Alcubierre default left INCOMPLETE rays' - assert any(v[10] == 0 for v in alc_vertices), \ - 'Alcubierre default produced no escaped ray' - assert sum_rhs(alc_vertices) > 0, \ - 'Alcubierre default pre-routed every ray; lookback misses feature' - assert image_payload(alc_out) + assert '|v_s| < 1' not in alc_help, 'help still claims |v_s| < 1' + escape1 = 1.0 + 20.0 / 1.0 # R = sigma = 1 - # The migrated coordinate-time coverage is the physical geometry budget - # margin*4*escape/(1-|v_s|) with escape = R + 20/sigma, and it is - # independent of the accepted-step count and of the chosen initial - # step. Two cheap maps with different resource overrides must keep the - # same lookback. - alc_escape = 1.0 + 20.0 / 1.0 # R + 20/sigma for R=sigma=1 - alc_sep = 1.0 - abs(0.3) - alc_expected_lookback = 1.25 * 4.0 * alc_escape / alc_sep - assert abs(alc_prov['max_lookback_time'] - alc_expected_lookback) \ - < 1e-12, alc_prov['max_lookback_time'] + def alc_camera(vs, ra_deg, dec_deg=0.0): + return ['--observer-position', '0', '0', '0', + '--observer-velocity', repr(vs), '0', '0', + '--look-ra-deg', repr(ra_deg), + '--look-dec-deg', repr(dec_deg)] - def alc_map_with(tag, *options): - m = tmp / f'alcubierre_{tag}.grlens' - run(alc, '--alcubierre-vs', 0.3, '--alcubierre-radius', 1, - '--alcubierre-sigma', 1, '--catalog', 'assets/sky_grid_5deg.csv', - '--width', 16, '--height', 8, '--fov-deg', 80, '--exposure', - 1e-3, '--coarse-cell-pixels', 8, '--refine-max-level', 0, - '--psf-relative-tail', 1e-4, '--allow-incomplete', - '--output', tmp / f'alcubierre_{tag}.{ext}', - '--lens-map-output', m, *options) - return map_provenance(m)[2] + def alc_map_prov(tag, vs, D=8.0, radius=1.0, extra=(), allow=True, + width=8, height=4, cell=4, camera=None, + catalog=True): + m = tmp / f'alc_{tag}.grlens' + args = ['--alcubierre-vs', repr(vs), '--alcubierre-radius', + repr(radius), '--alcubierre-sigma', '1', + '--dark-threshold', repr(D), '--width', str(width), + '--height', str(height), '--fov-deg', 80, '--exposure', + '1e-3', '--coarse-cell-pixels', str(cell), + '--refine-max-level', 0, '--psf-relative-tail', 1e-4] + if catalog: + args += ['--catalog', 'assets/sky_grid_5deg.csv'] + if camera is not None: + args += camera + if allow: + args.append('--allow-incomplete') + args += ['--output', str(tmp / f'alc_{tag}.{ext}'), + '--lens-map-output', str(m), *extra] + run(alc, *args) + return map_provenance(m)[2], m - steps_override = alc_map_with('steps_override', '--trace-max-steps', 8) - assert steps_override['initial_max_steps'] == 8 - assert abs(steps_override['max_lookback_time'] - alc_expected_lookback) \ - < 1e-12, steps_override - step_override = alc_map_with('step_override', '--ode-initial-step', 0.02) - assert abs(step_override['coordinate_time_step'] - 0.02) < 1e-15 - assert abs(step_override['max_lookback_time'] - alc_expected_lookback) \ - < 1e-12, step_override + # 1) Ordinary sub-luminal separation: B is exactly 5*escape/sep and is + # independent of the step count and of the initial step. + p_sub, _ = alc_map_prov('sub', 0.3, camera=alc_camera(0.3, 0.0)) + assert p_sub['integrator'] == 1, p_sub + assert math.isclose(p_sub['max_lookback_time'], + alcubierre_budget(escape1, 0.3, 8.0), rel_tol=1e-12) + assert p_sub['min_step'] <= p_sub['coordinate_time_step'] \ + <= p_sub['max_step'] + p_steps, _ = alc_map_prov('sub_steps', 0.3, + extra=('--trace-max-steps', '8'), + camera=alc_camera(0.3, 0.0)) + assert p_steps['initial_max_steps'] == 8 + assert math.isclose(p_steps['max_lookback_time'], + alcubierre_budget(escape1, 0.3, 8.0), rel_tol=1e-12) + p_step, _ = alc_map_prov('sub_step', 0.3, + extra=('--ode-initial-step', '0.02'), + camera=alc_camera(0.3, 0.0)) + assert abs(p_step['coordinate_time_step'] - 0.02) < 1e-15 + assert math.isclose(p_step['max_lookback_time'], + alcubierre_budget(escape1, 0.3, 8.0), rel_tol=1e-12) - # Extreme separation (v_s = 0.99999999) where the legacy fixed-step - # estimate far exceeds the cap. The DP path must accept an explicit - # tiny coordinate-time budget instead of being rejected at startup by - # the fixed-step guard, and it must actually exercise the quota path - # (UNRESOLVED rays or real RHS work), not pre-route everything to - # escapes. --allow-incomplete publishes the diagnostic frame. - extreme = ['--alcubierre-vs', 0.99999999, '--alcubierre-radius', 1, - '--alcubierre-sigma', 1, '--catalog', - 'assets/sky_grid_5deg.csv', '--width', 8, '--height', 4, - '--fov-deg', 80, '--exposure', 1e-3, - '--coarse-cell-pixels', 4, '--refine-max-level', 0, - '--psf-relative-tail', 1e-4, - '--observer-position', 0, 0, 0, - '--observer-velocity', 0.99999999, 0, 0, - '--look-ra-deg', 0, '--look-dec-deg', 0] - ext_map = tmp / 'alcubierre_extreme.grlens' + # 2) Near-luminal (vs = 1, separation 0): the exp(-D) floor makes the + # allowance finite, and it grows with the dark threshold D. + p8, _ = alc_map_prov('near8', 1.0, D=8.0, camera=alc_camera(1.0, 0.0)) + p12, _ = alc_map_prov('near12', 1.0, D=12.0, camera=alc_camera(1.0, 0.0)) + assert math.isclose(p8['max_lookback_time'], + alcubierre_budget(escape1, 1.0, 8.0), rel_tol=1e-12) + assert math.isclose(p12['max_lookback_time'], + alcubierre_budget(escape1, 1.0, 12.0), rel_tol=1e-12) + assert p12['max_lookback_time'] > p8['max_lookback_time'] + + # 3) Both sides of the threshold-derived vcut and exactly luminal + # values use either separation or the finite floor; a tiny film traces + # and resolves without INCOMPLETE outcomes. + vcut = 1.0 - math.exp(-8.0) + for vs in (0.999, math.nextafter(vcut, 0.0), + math.nextafter(vcut, 1.0), 0.9999, math.nextafter(1.0, 0.0), + math.nextafter(1.0, 2.0)): + pv, mv = alc_map_prov(f'vcut_{vs!r}', vs, camera=alc_camera(vs, 0.0)) + assert math.isclose(pv['max_lookback_time'], + alcubierre_budget(escape1, vs, 8.0), + rel_tol=1e-12), (vs, pv['max_lookback_time']) + vv, _ = map_vertices(mv) + assert 3 not in {v[10] for v in vv}, (vs, 'INCOMPLETE ray') + assert sum_rhs(vv) > 0, vs + + # 4) The default camera (0,0,15 for R=5) must work with a superluminal + # bubble: the generic radius-15 camera lies inside escape radius 25. + pd, md = alc_map_prov('defaultcam', 2.0, radius=5.0) + dv, _ = map_vertices(md) + assert 3 not in {v[10] for v in dv}, 'default camera left INCOMPLETE rays' + assert sum_rhs(dv) > 0 + + # 5) Real small images at the critical velocities: the direction along + # the bubble motion is the DARK direction (Pi_x = sign(vs)); the + # central vertex must be DARK and some edge ray must escape. The + # 16x8 film with a 4-pixel coarse cell places a vertex exactly at the + # image center. + for vs in (1.0, -1.0, 2.0, -2.0): + ra = 0.0 if vs < 0.0 else 180.0 + pimg, mimg = alc_map_prov(f'img_{vs!r}', vs, width=16, height=8, + cell=4, camera=alc_camera(vs, ra)) + verts, _ = map_vertices(mimg) + outcomes = {v[10] for v in verts} + assert 3 not in outcomes, (vs, 'INCOMPLETE outcome', outcomes) + assert 0 in outcomes and 1 in outcomes, (vs, outcomes) + assert sum_rhs(verts) > 0, vs + central = min(verts, key=lambda v: (v[0] - 8.0) ** 2 + + (v[1] - 4.0) ** 2) + assert central[10] == 1, \ + (vs, 'central vertex is not DARK', central[0], central[1], + central[10]) + + # 6) Extreme tiny-budget quota path: an explicit 4-step, 1-time + # allowance from an inside camera must stop as UNRESOLVED, never as + # a fabricated escape from a trace that took no steps. + extreme = (alc_camera(0.99999999, 0.0) + + ['--alcubierre-vs', '0.99999999', '--alcubierre-radius', + '1', '--alcubierre-sigma', '1', '--catalog', + 'assets/sky_grid_5deg.csv', '--width', '8', '--height', + '4', '--fov-deg', '80', '--exposure', '1e-3', + '--coarse-cell-pixels', '4', '--refine-max-level', 0, + '--psf-relative-tail', '1e-4']) + ext_map = tmp / 'alc_extreme.grlens' run(alc, *extreme, '--integrator', 'dp54', - '--trace-lookback-time', 1, '--trace-max-steps', 4, - '--max-total-steps', 4, '--retry-step-increment', 0, - '--max-total-lookback-time', 1, '--retry-lookback-increment', 0, - '--allow-incomplete', '--output', - tmp / f'alcubierre_extreme.{ext}', '--lens-map-output', ext_map) + '--trace-lookback-time', '1', '--trace-max-steps', '4', + '--max-total-steps', '4', '--retry-step-increment', '0', + '--max-total-lookback-time', '1', '--retry-lookback-increment', '0', + '--allow-incomplete', '--output', tmp / f'alc_extreme.{ext}', + '--lens-map-output', ext_map) _, _, ext_prov = map_provenance(ext_map) assert ext_prov['integrator'] == 1 assert abs(ext_prov['max_lookback_time'] - 1.0) < 1e-15, ext_prov assert ext_prov['initial_max_steps'] == 4, ext_prov ext_vertices, _ = map_vertices(ext_map) - assert any(v[10] == 2 for v in ext_vertices) or \ - sum_rhs(ext_vertices) > 0, \ - 'extreme Alcubierre case did not exercise the quota path' + assert all(v[10] == 2 for v in ext_vertices), \ + ('a 1-time/4-step trace fabricated a non-UNRESOLVED outcome', + [v[10] for v in ext_vertices]) - # The same parameters without explicit DP budgets keep the legacy - # fixed-step startup guard, which still rejects the estimated domain. - rk4_extreme = run( - alc, '--integrator', 'rk4', '--alcubierre-vs', 0.99999999, - '--alcubierre-radius', 1, '--alcubierre-sigma', 1, '--catalog', - 'assets/sky_grid_5deg.csv', '--width', 8, '--height', 4, - '--fov-deg', 80, '--exposure', 1e-3, '--coarse-cell-pixels', 4, - '--refine-max-level', 0, '--psf-relative-tail', 1e-4, - '--output', tmp / f'alcubierre_rk4_extreme.{ext}', ok=False) - assert rk4_extreme.returncode == 2, rk4_extreme.stderr - assert 'cap' in rk4_extreme.stderr, rk4_extreme.stderr + # 7) Large dark threshold (D=1000) with vs=1: exp(-1000) underflows to + # 0 and the separation is exactly 0, so the log fallback saturates + # the default allowance to DBL_MAX/4. A 4-step/4-total cap with + # retry disabled keeps the trace short while the provenance records + # the non-masked saturated default. An explicit --trace-lookback + # overrides it without masking the independent step allowance. + psat, _ = alc_map_prov('sat1000', 1.0, D=1000.0, + extra=('--trace-max-steps', '4', + '--max-total-steps', '4', + '--retry-step-increment', '0'), + camera=alc_camera(1.0, 0.0)) + assert math.isclose(psat['max_lookback_time'], + sys.float_info.max / 4.0, rel_tol=1e-12), \ + psat['max_lookback_time'] + assert psat['initial_max_steps'] == 4 + psat_ov, _ = alc_map_prov('sat1000_ov', 1.0, D=1000.0, + extra=('--trace-lookback-time', '1', + '--trace-max-steps', '4', + '--max-total-steps', '4', + '--retry-step-increment', '0'), + camera=alc_camera(1.0, 0.0)) + assert psat_ov['max_lookback_time'] == 1.0, psat_ov + assert psat_ov['initial_max_steps'] == 4, psat_ov + # Explicit time override must not mask the (default) step allowance. + pind, _ = alc_map_prov('indep', 1.0, D=8.0, + extra=('--trace-lookback-time', '1'), + camera=alc_camera(1.0, 0.0)) + assert pind['max_lookback_time'] == 1.0, pind + assert pind['initial_max_steps'] > 1, pind - print('alcubierre: default DP config validated and produced escapes; ' - 'extreme DP quota path ok, RK4 guard still rejects', flush=True) + # 8) RK4 guard: a super-luminal default (vs=2, small estimate) is no + # longer rejected; the guard triggers only when the *derived default* + # estimate exceeds the cap (D=12, vs=1), and an explicit + # --trace-max-steps is a user allowance that bypasses it. + run(alc, '--integrator', 'rk4', '--alcubierre-vs', '2', + '--alcubierre-radius', '1', '--alcubierre-sigma', '1', + '--catalog', 'assets/sky_grid_5deg.csv', '--width', '8', + '--height', '4', '--fov-deg', '80', '--exposure', '1e-3', + '--coarse-cell-pixels', '4', '--refine-max-level', 0, + '--psf-relative-tail', '1e-4', '--allow-incomplete', + '--output', tmp / f'alc_rk4_v2.{ext}') + rk4_guard = run( + alc, '--integrator', 'rk4', '--alcubierre-vs', '1', + '--alcubierre-radius', '1', '--alcubierre-sigma', '1', + '--dark-threshold', '12', '--catalog', 'assets/sky_grid_5deg.csv', + '--width', '8', '--height', '4', '--fov-deg', '80', + '--exposure', '1e-3', '--coarse-cell-pixels', '4', + '--refine-max-level', 0, '--psf-relative-tail', '1e-4', + '--output', tmp / f'alc_rk4_guard.{ext}', ok=False) + assert rk4_guard.returncode == 2, rk4_guard.stderr + assert 'cap' in rk4_guard.stderr, rk4_guard.stderr + run(alc, '--integrator', 'rk4', '--alcubierre-vs', '1', + '--alcubierre-radius', '1', '--alcubierre-sigma', '1', + '--dark-threshold', '12', '--trace-max-steps', '4', + '--catalog', 'assets/sky_grid_5deg.csv', '--width', '8', + '--height', '4', '--fov-deg', '80', '--exposure', '1e-3', + '--coarse-cell-pixels', '4', '--refine-max-level', 0, + '--psf-relative-tail', '1e-4', '--allow-incomplete', + '--output', tmp / f'alc_rk4_expl.{ext}') + + # 9) Impossible parameters stay clearly rejected. + for bad in (['--alcubierre-vs', 'nan'], ['--alcubierre-vs', 'inf'], + ['--alcubierre-radius', '0'], ['--alcubierre-sigma', '0']): + bad_out = tmp / f'alc_bad.{ext}' + rejected = run(alc, *bad, '--catalog', 'assets/sky_grid_5deg.csv', + '--width', '8', '--height', '4', '--fov-deg', '80', + '--exposure', '1e-3', '--coarse-cell-pixels', '4', + '--refine-max-level', 0, '--psf-relative-tail', '1e-4', + '--output', bad_out, ok=False) + assert rejected.returncode != 0, rejected.stderr + assert not bad_out.exists() + + print('alcubierre: budget/provenance, vcut scans, vcut images with ' + 'central DARK, default camera, large-D saturation, explicit ' + 'overrides and RK4 guard all passed', flush=True) diff --git a/tests/test_alcubierre.c b/tests/test_alcubierre.c index a63d4c6..df0032a 100644 --- a/tests/test_alcubierre.c +++ b/tests/test_alcubierre.c @@ -35,8 +35,24 @@ int main(void) { SpacetimeSource source = {0}; MetricData metric; CHECK(spacetime_create_alcubierre(&source, vs, radius, sigma) == 0); - CHECK(spacetime_create_alcubierre(&source, 1.0, radius, sigma) != 0); - CHECK(spacetime_create_alcubierre(&source, -1.5, radius, sigma) != 0); + /* Sub- and super-luminal velocities are both accepted. A successful + * constructor installs a context, so each acceptance uses its own temporary + * source that is destroyed immediately; the shared `source` above is never + * overwritten with a second live context. */ + { + SpacetimeSource luminal = {0}; + CHECK(spacetime_create_alcubierre(&luminal, 1.0, radius, sigma) == 0); + spacetime_destroy(&luminal); + } + { + SpacetimeSource superluminal = {0}; + CHECK(spacetime_create_alcubierre(&superluminal, -1.5, radius, sigma) == 0); + spacetime_destroy(&superluminal); + } + /* Non-finite velocities stay rejected. */ + CHECK(spacetime_create_alcubierre(&source, NAN, radius, sigma) != 0); + CHECK(spacetime_create_alcubierre(&source, INFINITY, radius, sigma) != 0); + CHECK(spacetime_create_alcubierre(&source, -INFINITY, radius, sigma) != 0); CHECK(spacetime_create_alcubierre(&source, vs, 0.0, sigma) != 0); CHECK(spacetime_create_alcubierre(&source, vs, radius, 0.0) != 0); /* A derived escape radius that overflows or does not exceed R is rejected. */ @@ -322,6 +338,99 @@ int main(void) { } } + /* Production DP54 axial superluminal check. A comoving bubble-center camera + * at R = 1, sigma = 1 (escape radius 21) sees Pi_x = +-1 along the bubble + * axis. The shared camera-relative dark policy must fire at threshold 8 for + * the direction along the bubble motion and the opposite direction must + * escape. The vs = +-2 constants were computed independently with mpmath; + * this test has no dependency on any local experiment fixture. */ + { + const double R1 = 1.0, sig1 = 1.0; + static const double vlist[] = {1.0, -1.0, 2.0, -2.0, 0.9999}; + for (size_t k = 0; k < sizeof vlist / sizeof vlist[0]; ++k) { + const double v = vlist[k]; + const double sgn = v > 0.0 ? 1.0 : -1.0; + SpacetimeSource fast = {0}; + CHECK(spacetime_create_alcubierre(&fast, v, R1, sig1) == 0); + ObserverCamera cam = {.coordinate_time = 0.0, + .position = {0.0, 0.0, 0.0}, + .velocity = {v, 0.0, 0.0}, + .look_ra_deg = 0.0, + .look_dec_deg = 0.0, + .roll_deg = 0.0}; + MetricData m; + CHECK(eval(&fast, 0.0, cam.position, &m) == 0); + ObserverState o; + CHECK(observer_from_coordinate_camera(&m, &cam, &o, NULL) == + OBSERVER_BUILD_OK); + /* A coordinate-static center camera is timelike only for |v| < 1. */ + ObserverCamera stat = cam; + stat.velocity[0] = stat.velocity[1] = stat.velocity[2] = 0.0; + ObserverState so; + const int static_ok = observer_from_coordinate_camera(&m, &stat, &so, + NULL) == + OBSERVER_BUILD_OK; + CHECK(static_ok == (fabs(v) < 1.0)); + /* An outer static camera in the flat exterior is always legal. */ + { + const double pos[3] = {26.0, 0.0, 0.0}; + MetricData om; + ObserverCamera oc = {.coordinate_time = 0.0, + .position = {26.0, 0.0, 0.0}, + .look_ra_deg = 0.0, + .look_dec_deg = 0.0, + .roll_deg = 0.0}; + CHECK(eval(&fast, 0.0, pos, &om) == 0); + ObserverState oo; + CHECK(observer_from_coordinate_camera(&om, &oc, &oo, NULL) == + OBSERVER_BUILD_OK); + } + const GeodesicTraceConfig trace = { + .coordinate_time_step = 0.05, + .max_steps = 100000u, + .threshold = {.kind = THRESHOLD_LOG_ENERGY_GROWTH, + .value = 8.0, + .policy_version = 1}, + .stepper = GEODESIC_STEPPER_DP54, + .atol_x = 1e-9, + .atol_Pi = 1e-9, + .atol_L = 1e-9, + .rtol = 1e-9, + .min_step = 1e-12, + .max_step = 0.4, + .consecutive_rejection_limit = 32, + .max_lookback_time = 30000.0}; + const double n_dark[3] = {-sgn, 0.0, 0.0}; + const double n_esc[3] = {sgn, 0.0, 0.0}; + GeodesicRayState idark, iesc; + CHECK(geodesic_initialize_past_ray_metric(&m, &o, n_dark, &idark) == 0); + CHECK(geodesic_initialize_past_ray_metric(&m, &o, n_esc, &iesc) == 0); + printf("alcubierre vs=%.6g dark Pi_x=%.17g escape Pi_x=%.17g\n", v, + idark.Pi[0], iesc.Pi[0]); + CHECK(sgn * idark.Pi[0] > 0.999 && sgn * idark.Pi[0] < 1.000000001); + CHECK(sgn * iesc.Pi[0] < -0.999 && sgn * iesc.Pi[0] > -1.000000001); + const RayEndpoint dark = geodesic_trace_past(&fast, &o, n_dark, &trace); + CHECK(dark.outcome == RAY_OUTCOME_DARK); + CHECK(isfinite(dark.stop_coordinate_time)); + CHECK(fabs(dark.threshold_value - 8.0) < 1e-6); + CHECK(fabs((dark.final_log_alpha_p0 - dark.final_log_alpha_p0_0) - 8.0) < + 1e-6); + if (fabs(v) == 2.0) { + const double q = dark.final_x[0] - v * dark.stop_coordinate_time; + CHECK(fabs(fabs(q) - 1.2181434100155241) < 1e-6); + CHECK(fabs(dark.stop_coordinate_time + 7.09915163394274) < 1e-6); + } + const RayEndpoint esc = geodesic_trace_past(&fast, &o, n_esc, &trace); + CHECK(esc.outcome == RAY_OUTCOME_ESCAPED); + CHECK(isfinite(esc.frequency_ratio) && esc.frequency_ratio > 0.0); + if (fabs(v) == 2.0) + CHECK(fabs(esc.frequency_ratio - 3.0) < 1e-6); + for (int i = 0; i < 3; ++i) + CHECK(isfinite(esc.n_infinity[i])); + spacetime_destroy(&fast); + } + } + spacetime_destroy(&source); puts("alcubierre regression passed"); return 0; diff --git a/usage.md b/usage.md index a6d5fff..f6777d6 100644 --- a/usage.md +++ b/usage.md @@ -116,32 +116,39 @@ with the bubble center following the constant-velocity worldline $$f(r) = \frac{\tanh(\sigma(r+R)) - \tanh(\sigma(r-R))}{2\tanh(\sigma R)},\qquad r_s = \sqrt{(x-x_s)^2 + y^2 + z^2}.$$ -The bubble therefore propagates through the coordinates, and the metric is -time-dependent: the renderer evaluates `f(r_s)` and its spatial derivatives at -each coordinate time, while the extrinsic curvature supplies the required -`d_t gamma` information to the 3+1 null-ray equations. The exotic matter that -would source the bubble is treated as optically transparent and there is no -horizon, so in practice rays are active or escaped: the bubble has no causal -boundary at which `L - L0` can diverge, and the shared dark policy is not -expected to trigger. This is why the backend requires a sub-luminal -`|v_s| < 1`; at or above `1` the metric develops an ergoregion/event horizon and -static observers cease to exist, which is outside the current scope. +The renderer evaluates the moving bubble's time-dependent metric along each +ray. The exotic matter sourcing the bubble is treated as optically transparent. +The shared dark policy terminates rays when `L - L0 >= T` (`T = 8` by default, +set with `--dark-threshold`); this finite threshold can also be reached at +sub-luminal bubble velocities. | Option | Meaning / default | | --- | --- | -| `--alcubierre-vs V` | Constant shift parameter, `|V| < 1` (default 0.5) | +| `--alcubierre-vs V` | Constant bubble velocity `v_s` (any finite value, default 0.5) | | `--alcubierre-radius R` | Bubble radius `R > 0` (default 5) | | `--alcubierre-sigma S` | Wall sharpness `S > 0` (default 1) | -`f` decays to zero past `r_s = R` over a transition width `~1/sigma`, so the -finite escape sphere is bubble-centered with radius `R + 20/sigma` and needs no -CLI option; it follows the moving bubble, so rays terminate only once the local -metric is flat to below double precision. The single-frame camera default is -`(0,0,15)` at `t = 0`, when the bubble is still at the origin; it must lie -inside the escape sphere, or the observer build fails with an explicit error. -The per-ray step budget scales with the escape radius and `1/(1-|v_s|)`, so -near-luminal `v_s` still lets grazing rays escape; combinations whose -worst-case budget would exceed the internal cap are rejected at startup. +A camera must be timelike with an orthonormal tetrad. Its coordinate velocity +`V = dx/dt` must satisfy `|V - v_s f e_x| < 1`. A static camera requires +`|v_s f| < 1`; at the bubble center, the comoving velocity `(v_s, 0, 0)` is +timelike even for super-luminal bubbles. Cameras may lie inside or outside the +escape sphere; exterior rays are routed to their first entry or to infinity. + +The escape sphere follows the bubble with radius `R + 20/sigma`; escaping rays +continue through a Minkowski exterior. The single-frame camera defaults to +`(0,0,15)` at `t = 0`. The finite-radius truncation leaves a residual shift of +order `|v_s| e^{-40}`; accuracy at extremely large velocities is not guaranteed. + +Per-ray coordinate-time coverage is a **resource allowance** + +$$B = \frac{5\,R_\text{escape}}{\max\!\big(|1-|v_s||,\; e^{-T}\big)},\qquad +R_\text{escape} = R + \frac{20}{\sigma},$$ + +with `T = --dark-threshold`, clamped to `[DBL_MIN, DBL_MAX/4]`. This is not a +completion guarantee: quota exhaustion returns `UNRESOLVED/BUDGET_EXHAUSTED`. +`--trace-lookback-time` overrides `B` independently of `--trace-max-steps`. +RK4 rejects a derived default step estimate above its internal cap; an explicit +`--trace-max-steps` bypasses that check. Lensing and frequency shifts come from the bubble wall. The configuration is invariant under the isometry `(t, x) -> (t + T, x + v_s T)`, so observers @@ -152,7 +159,7 @@ example: make -j PSF_BACKEND=cpu SPACETIME=alcubierre backend ./build/Release/alcubierre_sky --catalog assets/sky_grid_5deg.csv \ --observer-radius 15 --look-ra-deg 90 --look-dec-deg -90 \ - --alcubierre-vs 0.5 --alcubierre-radius 5 --alcubierre-sigma 1 \ + --alcubierre-vs 1.5 --alcubierre-radius 1 --alcubierre-sigma 1 \ --width 640 --height 360 --fov-deg 60 --exposure 1 \ --coarse-cell-pixels 16 --refine-max-level 2 --psf-direct \ --output output/imgs/alcubierre_wall.png