diff --git a/.gitignore b/.gitignore index 92f220a..9bec68d 100644 --- a/.gitignore +++ b/.gitignore @@ -28,3 +28,5 @@ tests/coverage_html/ # python tests/fetch_data.py # once after cloning to populate this directory. Unit tests don't need it. tests/sxs_cache/ + +dist/ \ No newline at end of file diff --git a/docs/source/index.rst b/docs/source/index.rst index 0f5b20d..27ba066 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -17,6 +17,7 @@ Features of gwBOB - Generate any flavor of BOB (psi4, news, strain, mass and current quadrupoles) using various initial conditions! - Easy calling of SXS and CCE waveforms as well as user provided NR waveforms! +- Standalone mode: build a BOB waveform from just :math:`(m_f, \chi_f, l, m)` — no NR data required! - Easy comparisons of BOB to NR waveforms as well as a sum of overtones! - "Simple" parameter estimation using BOB! - Open source, documented, and actively developed! diff --git a/docs/source/quickstart.rst b/docs/source/quickstart.rst index f7cdaa9..838662f 100644 --- a/docs/source/quickstart.rst +++ b/docs/source/quickstart.rst @@ -124,3 +124,37 @@ Integrating the news to obtain the strain is a complex problem separate from the .. image:: images/BOB_strain_0305.png + +Standalone mode (no NR data required) +------------------------------------- + +Sometimes you just want a BOB waveform for a remnant whose final mass and spin you already know. ``initialize_standalone`` does exactly that. + +.. code-block:: python + + BOB = BOB_utils.BOB() + BOB.initialize_standalone(mf=0.95, chif=0.7) + BOB.what_should_BOB_create = "news" + t, y = BOB.construct_BOB() + +The peak time defaults to :math:`t_p = 0`, the peak amplitude defaults to :math:`A_p = 1` , the QNM parameters are computed from Kerr via the ``qnm`` package, and :math:`\Omega_0` is set to the mode-appropriate fit value (``Omega_0_fit_psi4``, ``Omega_0_fit_news``, or ``Omega_0_fit_strain``) at the moment you set ``what_should_BOB_create``. Switching modes refits :math:`\Omega_0` automatically. + +You can also pass any of these explicitly: + +.. code-block:: python + + BOB.initialize_standalone( + mf=0.95, chif=0.7, l=2, m=2, + tp=0.0, Ap=1.0, # peak time and amplitude + start_before_tpeak=-75.0, + end_after_tpeak=75.0, + resample_dt=0.1, # output grid spacing (>= 0.1) + Omega_0=0.123, # explicit Omega_0; sticky across mode switches + t0=-10.0, # finite-t0 build, relative to tp + w_r=0.5, tau=12.0, # non-Kerr QNM (e.g., for "Beyond Kerr" experiments) + ) + +Because no NR data is loaded, several capabilities are unavailable in standalone mode and raise a clear error if you try them: ``optimize_Omega0`` and the other ``optimize_*`` flags, the ``fit_*`` drivers, the ``construct_BOB_*_quadrupole_naturally`` helpers, the per-mode ``get_psi4_data`` / ``get_news_data`` / ``get_strain_data`` accessors, and the ``strain_using_news`` / ``strain_using_psi4`` / ``mass_quadrupole_*`` / ``current_quadrupole_*`` modes (because these need NR timeseries for the integration constants or for the :math:`(l, \pm m)` combination). + +Supported modes in standalone mode are ``"psi4"``, ``"news"``, and ``"strain"``. The finite-:math:`t_0` build is available — pass ``t0=`` at init time as shown above, or use ``BOB.set_initial_time = -10`` after setting ``what_should_BOB_create``. ``BOB.NR_based_on_BOB_ts`` is ``None`` after ``construct_BOB`` because there is no NR data to align against; use the returned ``(t, y)`` arrays directly. + diff --git a/gwBOB/BOB_utils.py b/gwBOB/BOB_utils.py index 54dc3ca..76e9acb 100644 --- a/gwBOB/BOB_utils.py +++ b/gwBOB/BOB_utils.py @@ -51,6 +51,16 @@ "current_quadrupole_with_psi4": ("current_quadrupole_with_psi4", "psi4", "current"), } +# Claude Code: Dispatch table for standalone (no-NR) modes. See DESIGN_standalone_init.md. +# Trimmed subset of _SIMPLE_MODES — only the three modes that don't need NR data: +# strain_using_* needs the NR series-integration constants and quadrupole modes +# need NR (l, ±m) timeseries. Tuple is (canonical_name, Omega_0 fit fn). +_STANDALONE_MODES = { + "psi4": ("psi4", gen_utils.Omega_0_fit_psi4), + "news": ("news", gen_utils.Omega_0_fit_news), + "strain": ("strain", gen_utils.Omega_0_fit_strain), +} + class BOB: ''' @@ -151,7 +161,9 @@ def what_should_BOB_create(self, value): val = value.lower() self._ensure_qnm_data_ready() - if val in _SIMPLE_MODES: + if self._runtime.is_standalone: + self._apply_standalone_mode(val) + elif val in _SIMPLE_MODES: self._apply_simple_mode(val) elif val in _QUADRUPOLE_MODES: self._apply_quadrupole_mode(val) @@ -219,6 +231,36 @@ def _apply_quadrupole_mode(self, val): else: self._runtime.tp = tp + def _apply_standalone_mode(self, val): + '''Internal: handle the three psi4/news/strain modes after standalone init. + + Two deliberate differences from ``_apply_simple_mode``: + + 1. ``tp`` and ``Ap`` are NOT recomputed — they were captured at + ``initialize_standalone`` time from user input, not derived from NR + data. The trailing time-grid recompute in the setter still uses + ``self._runtime.tp`` correctly because it was set at init. + + 2. The ``Omega_0`` fit is gated on ``Omega0_user_override`` so an + explicit user-supplied ``Omega_0`` survives mode switches. + + Claude Code: See DESIGN_standalone_init.md "Behavioral changes elsewhere". + ''' + if val not in _STANDALONE_MODES: + raise ValueError( + f"In standalone mode, what_should_BOB_create must be one of " + f"{sorted(_STANDALONE_MODES)}. The mode {val!r} is not supported " + "because it requires NR data. Re-initialize with " + "initialize_with_sxs_data, initialize_with_cce_data, or " + "initialize_with_NR_mode if you need this capability." + ) + canonical, omega_fn = _STANDALONE_MODES[val] + self._wf_config.what_to_create = canonical + self._runtime.data = None # explicit: no NR data exists in standalone mode + # _runtime.tp and _runtime.Ap were set by initialize_standalone — leave them. + if not self._runtime.Omega0_user_override: + self._remnant.Omega_0 = omega_fn(self._remnant.mf, self._remnant.chif_with_sign) + @property def set_initial_time(self): '''Return the configured initial time t0 (in code units relative to peak).''' @@ -232,6 +274,14 @@ def set_initial_time(self,value): float, the initial time is set to the value and the frequency is set using the data specified by "what_should_BOB_create". + Standalone-mode behavior (Claude Code: see DESIGN_standalone_init.md): + when ``self._runtime.is_standalone`` is True, the NR-frequency lookup + is skipped — there is no NR data to read. ``t0`` and ``minf_t0`` are + set as usual; ``Omega_0`` retains whatever value it had (the + mode-appropriate fit applied at ``what_should_BOB_create`` time, or + the user override). The user can set ``Omega_0`` manually if they + want a different value for the finite-t0 build. + args: value (tuple or float): Initial time and whether to set the frequency using the strain data ''' @@ -245,6 +295,19 @@ def set_initial_time(self,value): set_freq_using_strain_data = False self._wf_config.minf_t0 = False + if self._runtime.is_standalone: + # No NR data — skip the Omega_0 refit entirely. t0 is interpreted + # in the same coordinate as the NR-init path: relative to tp. + if set_freq_using_strain_data: + logger.warning( + "set_initial_time: tuple form (use strain frequency) is " + "ignored in standalone mode — no NR strain data exists. " + "Omega_0 left at its current value." + ) + self._wf_config.t0 = self._runtime.tp + value + self._runtime.t0_tp_tau = (self._wf_config.t0 - self._runtime.tp)/self._remnant.tau + return + if(set_freq_using_strain_data): freq = gen_utils.get_frequency(self._data.strain_data) else: @@ -320,6 +383,8 @@ def optimize_t0_and_Omega0(self,value): args: value (bool): Optimize initial time and frequency ''' + if value: + self._require_NR("optimize_t0_and_Omega0") self._wf_config.minf_t0 = False self._fit_config.optimize_t0_and_Omega0 = value @@ -337,6 +402,8 @@ def optimize_t0(self,value): args: value (bool): Optimize initial time ''' + if value: + self._require_NR("optimize_t0") self._wf_config.minf_t0 = False self._fit_config.optimize_t0 = value @@ -498,7 +565,10 @@ def what_is_BOB_building(self, v): self._wf_config.what_is_BOB_building = v @property def optimize_Omega0(self): return self._fit_config.optimize_Omega0 @optimize_Omega0.setter - def optimize_Omega0(self, v): self._fit_config.optimize_Omega0 = v + def optimize_Omega0(self, v): + if v: + self._require_NR("optimize_Omega0") + self._fit_config.optimize_Omega0 = v @property def start_fit_before_tpeak(self): return self._fit_config.start_fit_before_tpeak @start_fit_before_tpeak.setter @@ -670,6 +740,7 @@ def fit_omega(self,x,Omega_0): np.ndarray: BOB frequency evaluated at ``self.t[start:end]`` for the configured fit window. ''' + self._require_NR("fit_omega") #this function can be called if X_using_Y. self.Omega_0 = Omega_0 if('psi4' in self._wf_config.what_to_create): @@ -705,6 +776,7 @@ def fit_t0_and_omega(self,x,t0,Omega_0): evaluation failure, returns a flat array of 1e10 to signal a bad residual to the optimizer. ''' + self._require_NR("fit_t0_and_omega") #this function can be called if X_using_Y. self.Omega_0 = Omega_0 self.t0 = t0 @@ -745,6 +817,7 @@ def residual_t0_and_omega(self,p,t_freq,y_freq): float: Sum of squared residuals between the BOB-predicted and NR frequencies over the configured fit window. ''' + self._require_NR("residual_t0_and_omega") #freq = gen_utils.get_frequency(self.data) freq = kuibit_ts(t_freq,y_freq) t0,Omega_0 = p @@ -789,8 +862,9 @@ def fit_t0_only(self,t00,freq_data): float: Sum of squared residuals between BOB and NR frequencies on the configured fit window. ''' + self._require_NR("fit_t0_only") #freq data passed in is big Omega, where w = m*Omega - self.t0 = t00[0] + self.t0 = t00[0] self.t0_tp_tau = (self.t0 - self.tp)/self.tau self.Omega_0 = freq_data.y[gen_utils.find_nearest_index(freq_data.t,self.t0)] #freq data is already big Omega start_index = gen_utils.find_nearest_index(self.t,self.tp+self.start_fit_before_tpeak) @@ -828,6 +902,7 @@ def fit_Omega0(self): Raises: ValueError: If ``minf_t0`` is False (use ``fit_t0`` for finite t0). ''' + self._require_NR("fit_Omega0") if(self.minf_t0 is False): raise ValueError("You are setup for a finite t0 right now. Omega_0 fitting is only defined for t0 = infinity.") if(self._wf_config.end_after_tpeak 0.01 or abs(chif_arr[1]) > 0.01: + raise ValueError("Final spin has non-zero x or y component; not supported") + sign = float(np.sign(chif_arr[2])) or 1.0 + chif_mag = float(np.linalg.norm(chif_arr)) + else: + raise ValueError("chif must be a scalar or length-3 array") + + self.mf = mf + self.chif = chif_mag + self.chif_with_sign = sign * chif_mag + self.l = l + self.m = m + + self.Omega_ISCO = np.abs(gen_utils.get_Omega_isco(self.mf, self.chif_with_sign)) + # Claude Code: see DESIGN_standalone_init.md "Omega_0 defaulting". + # If user supplied Omega_0, store and mark override; else install + # Omega_ISCO as a placeholder (mirrors initialize_with_NR_mode) which + # _apply_standalone_mode overwrites with the mode-appropriate fit when + # what_should_BOB_create is set. + if Omega_0 is None: + self.Omega_0 = self.Omega_ISCO + self._runtime.Omega0_user_override = False + else: + self.Omega_0 = float(Omega_0) + self._runtime.Omega0_user_override = True + + if w_r < 0 or tau < 0: + logger.info("Calculating Kerr QNM parameters from provided mf and chif") + w_r_kerr, tau_kerr = gen_utils.get_qnm(self.mf, self.chif, self.l, np.abs(self.m), n=0, sign=sign) + self.w_r = np.abs(w_r_kerr) + self.tau = np.abs(tau_kerr) + else: + logger.info("Using user provided w_r and tau!") + self.w_r = np.abs(w_r) + self.tau = np.abs(tau) + self.Omega_QNM = self.w_r / np.abs(self.m) + + # Standalone-specific runtime state. Setting is_standalone last so any + # earlier failure leaves the BOB in a non-standalone (and thus + # consistent default) state. + self._runtime.tp = tp + self._runtime.Ap = Ap + self._data.resample_dt = resample_dt + self._wf_config.start_before_tpeak = start_before_tpeak + self._wf_config.end_after_tpeak = end_after_tpeak + if t0 is not None: + # Flip to the finite-t0 build. t0 is interpreted relative to tp, + # matching set_initial_time's coordinate convention. t0_tp_tau is + # also pre-computed here so the user can go straight to + # construct_BOB() after setting what_should_BOB_create. + self._wf_config.t0 = self._runtime.tp + float(t0) + self._wf_config.minf_t0 = False + self._runtime.t0_tp_tau = (self._wf_config.t0 - self._runtime.tp) / self.tau + self._runtime.is_standalone = True + + if verbose: + logger.debug("Standalone init: (l, m) = (%s, %s)", self.l, self.m) + logger.debug("Omega_ISCO = %s", self.Omega_ISCO) + logger.debug("Omega_QNM = %s", self.Omega_QNM) + logger.debug("tau = %s", self.tau) + logger.debug("tp = %s, Ap = %s", tp, Ap) + logger.debug("Omega_0 = %s (user override = %s)", + self.Omega_0, self._runtime.Omega0_user_override) + def _resolve_lm(self, kwargs): '''Return ``(l, m)`` from ``kwargs``, falling back to ``self.l`` / ``self.m``.''' l = kwargs.get('l', self.l) @@ -1718,6 +1943,34 @@ def _no_full_data_error(self, kind, l, m): "(memory-heavy — see MEMORY.md) to access non-default (l, m) modes." ) + def _require_NR(self, capability): + ''' + Raise ``RuntimeError`` if this BOB instance was initialized via the + standalone path (no NR data loaded). Called at the top of methods and + setters that read NR timeseries — fit drivers, optimize-flag setters, + per-mode data getters, and the quadrupole construction helpers. + + In non-standalone mode this is a no-op. + + Claude Code: See DESIGN_standalone_init.md "Disabled capabilities" for + the audit list and the rationale for centralizing the gate here. + + args: + capability (str): Human-readable name of the operation being + guarded; included in the error message so the user can locate + the call site. + + Raises: + RuntimeError: If ``self._runtime.is_standalone`` is True. + ''' + if self._runtime.is_standalone: + raise RuntimeError( + f"{capability} requires NR data; not available after " + f"initialize_standalone(). Use initialize_with_sxs_data, " + f"initialize_with_cce_data, or initialize_with_NR_mode if you " + f"need this capability." + ) + def get_psi4_data(self,**kwargs): ''' Return the NR psi4 timeseries for an arbitrary (l, m) mode. @@ -1740,6 +1993,7 @@ def get_psi4_data(self,**kwargs): returns: tuple of (np.ndarray, np.ndarray): (t, y) of the requested psi4 mode. ''' + self._require_NR("get_psi4_data") l, m = self._resolve_lm(kwargs) if self.full_psi4_data is None: if (l, m) == (self.l, self.m): return self.psi4_data.t, self.psi4_data.y @@ -1763,6 +2017,7 @@ def get_news_data(self,**kwargs): returns: tuple of (np.ndarray, np.ndarray): (t, y) of the requested news mode. ''' + self._require_NR("get_news_data") l, m = self._resolve_lm(kwargs) if self.full_strain_data is None: if (l, m) == (self.l, self.m): return self.news_data.t, self.news_data.y @@ -1785,6 +2040,7 @@ def get_strain_data(self,**kwargs): returns: tuple of (np.ndarray, np.ndarray): (t, y) of the requested strain mode. ''' + self._require_NR("get_strain_data") l, m = self._resolve_lm(kwargs) if self.full_strain_data is None: if (l, m) == (self.l, self.m): return self.strain_data.t, self.strain_data.y diff --git a/gwBOB/_state.py b/gwBOB/_state.py index 108a4a1..2f69096 100644 --- a/gwBOB/_state.py +++ b/gwBOB/_state.py @@ -102,6 +102,13 @@ class RuntimeState: t_tp_tau: Optional[np.ndarray] = None t0_tp_tau: Optional[np.ndarray] = None NR_based_on_BOB_ts: Any = None + # Claude Code: See DESIGN_standalone_init.md. Set by initialize_standalone; + # default False preserves all existing NR-init behavior. + is_standalone: bool = False + # Claude Code: See DESIGN_standalone_init.md. True when the user passed an + # explicit Omega_0 to initialize_standalone, signalling that + # _apply_standalone_mode must NOT overwrite it with the mode-appropriate fit. + Omega0_user_override: bool = False @dataclass(slots=True) diff --git a/tests/unit/test_BOB_standalone.py b/tests/unit/test_BOB_standalone.py new file mode 100644 index 0000000..20d98fc --- /dev/null +++ b/tests/unit/test_BOB_standalone.py @@ -0,0 +1,341 @@ +"""Unit tests for ``BOB.initialize_standalone`` (the no-NR-data path). + +Claude Code: See DESIGN_standalone_init.md for the architectural spec. + +Most tests here run without any NR data — that's the entire point of the +standalone path. The single exception is ``test_parity_with_NR_init``, +which uses the existing ``initial_sxs_bob_2325`` fixture from conftest.py +to confirm the standalone analytic minf-t0 build is identical to the +NR-init build for matching ``(mf, chif, tp, Ap, Omega_QNM, tau)``. +That test skips cleanly when the SXS cache is not populated. +""" + +from __future__ import annotations + +import copy + +import numpy as np +import pytest + +from gwBOB import gen_utils +from gwBOB.BOB_utils import BOB + + +# Default test parameters — light remnant, moderate prograde spin. +MF = 0.95 +CHIF = 0.7 + + +# --------------------------------------------------------------------------- +# Test 1: Happy path — psi4, news, strain. +# --------------------------------------------------------------------------- + +@pytest.mark.parametrize("mode", ["psi4", "news", "strain"]) +def test_happy_path(mode): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF, l=2, m=2) + bob.what_should_BOB_create = mode + t, y = bob.construct_BOB() + + assert t.shape == y.shape + assert y.dtype == np.complex128 + assert np.all(np.isfinite(y)) + + # Peak |y| must lie at t = 0 (the default tp = 0.0) within one grid spacing. + peak_t = t[np.argmax(np.abs(y))] + assert abs(peak_t) <= bob._data.resample_dt, ( + f"peak at t={peak_t}, expected ~0 within {bob._data.resample_dt}" + ) + + assert np.isfinite(bob.Omega_QNM) and bob.Omega_QNM > 0 + assert np.isfinite(bob.tau) and bob.tau > 0 + + +# --------------------------------------------------------------------------- +# Test 2: Convention parity with the NR-init path. +# --------------------------------------------------------------------------- + +def test_parity_with_NR_init(initial_sxs_bob_2325): + """A standalone BOB built with the NR-init BOB's tp/Ap/Omega_QNM/tau/Omega_0 + should produce the bit-identical analytic minf-t0 waveform. Proves no + numerical fork between the two code paths for matching inputs. + + Note: the NR path applies a phase alignment in ``construct_BOB``; the + standalone path skips it. We compare ``construct_BOB_minf_t0`` directly + to bypass that asymmetry — that's the analytic build proper. + """ + sxs_bob = copy.deepcopy(initial_sxs_bob_2325) + sxs_bob.what_should_BOB_create = "news" + sxs_ts = sxs_bob.construct_BOB_minf_t0(N=2) + + # Build a standalone BOB with the same physics constants. Pass Omega_0 + # explicitly so it matches sxs_bob.Omega_0 exactly (the standalone-mode + # default would re-fit and could differ at FP precision). + standalone = BOB() + standalone.initialize_standalone( + mf=sxs_bob.mf, + chif=sxs_bob.chif_with_sign, # signed scalar accepted + l=sxs_bob.l, + m=sxs_bob.m, + tp=sxs_bob.tp, + Ap=sxs_bob.Ap, + start_before_tpeak=sxs_bob._wf_config.start_before_tpeak, + end_after_tpeak=sxs_bob._wf_config.end_after_tpeak, + resample_dt=sxs_bob._data.resample_dt, + Omega_0=sxs_bob.Omega_0, + w_r=sxs_bob.w_r, + tau=sxs_bob.tau, + ) + standalone.what_should_BOB_create = "news" + standalone_ts = standalone.construct_BOB_minf_t0(N=2) + + np.testing.assert_allclose(sxs_ts.t, standalone_ts.t, rtol=1e-12, atol=0) + np.testing.assert_allclose(sxs_ts.y, standalone_ts.y, rtol=1e-12, atol=0) + + +# --------------------------------------------------------------------------- +# Test 3: chif shape acceptance (scalar, +z vector, -z vector). +# --------------------------------------------------------------------------- + +@pytest.mark.parametrize( + "chif_input, expected_signed", + [ + (0.7, 0.7), + (np.array([0.0, 0.0, 0.7]), 0.7), + (np.array([0.0, 0.0, -0.7]), -0.7), + ], +) +def test_chif_shape_acceptance(chif_input, expected_signed): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=chif_input) + # chif (unsigned magnitude) is always positive. + assert bob.chif == pytest.approx(0.7) + # chif_with_sign carries the sign of the z-component. + assert bob.chif_with_sign == pytest.approx(expected_signed) + + +# --------------------------------------------------------------------------- +# Test 4: Disabled-capability errors. Each must raise with a message +# mentioning "standalone" so the user can identify the constraint. +# --------------------------------------------------------------------------- + +def _fresh_standalone_bob(): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF) + bob.what_should_BOB_create = "news" + return bob + + +def test_disabled_optimize_Omega0(): + bob = _fresh_standalone_bob() + with pytest.raises(RuntimeError, match=r"(?i)standalone"): + bob.optimize_Omega0 = True + + +def test_disabled_strain_using_news(): + bob = _fresh_standalone_bob() + with pytest.raises(ValueError, match=r"(?i)standalone"): + bob.what_should_BOB_create = "strain_using_news" + + +def test_disabled_quadrupole_mode(): + bob = _fresh_standalone_bob() + with pytest.raises(ValueError, match=r"(?i)standalone"): + bob.what_should_BOB_create = "mass_quadrupole_with_news" + + +def test_disabled_get_psi4_data(): + bob = _fresh_standalone_bob() + with pytest.raises(RuntimeError, match=r"(?i)standalone"): + bob.get_psi4_data() + + +def test_disabled_fit_Omega0_direct_call(): + bob = _fresh_standalone_bob() + with pytest.raises(RuntimeError, match=r"(?i)standalone"): + bob.fit_Omega0() + + +def test_disabled_construct_BOB_current_quadrupole_naturally(): + bob = _fresh_standalone_bob() + with pytest.raises(RuntimeError, match=r"(?i)standalone"): + bob.construct_BOB_current_quadrupole_naturally() + + +# Walk every entry on the S3 audit list to confirm the _require_NR guard is +# wired in. Methods/properties that take no args are called with (); ones +# that need args get a minimal call shape that triggers the guard before +# the method does any real work. Some entry points are exercised via the +# named tests above; this parametrized test catches any guard that nobody +# else exercises so the audit list is locked in. +@pytest.mark.parametrize( + "operation", + [ + # Optimize-flag setters (the test_disabled_optimize_Omega0 case + # already covers Omega0; these check the other two). + lambda b: setattr(b, "optimize_t0", True), + lambda b: setattr(b, "optimize_t0_and_Omega0", True), + # Fit drivers and callbacks. + lambda b: b.fit_t0(), + lambda b: b.fit_t0_and_Omega0(), + lambda b: b.fit_omega(None, 0.1), + lambda b: b.fit_t0_and_omega(None, -10, 0.1), + lambda b: b.residual_t0_and_omega((-10, 0.1), None, None), + lambda b: b.fit_t0_only([-10], None), + # Quadrupole construction helpers. + lambda b: b.construct_NR_mass_and_current_quadrupole("news"), + lambda b: b.construct_BOB_mass_quadrupole_naturally(), + # Per-mode data getters. + lambda b: b.get_news_data(), + lambda b: b.get_strain_data(), + ], +) +def test_audit_list_complete(operation): + bob = _fresh_standalone_bob() + with pytest.raises(RuntimeError, match=r"(?i)standalone"): + operation(bob) + + +# --------------------------------------------------------------------------- +# Test 5: m=0 rejection. +# --------------------------------------------------------------------------- + +def test_m0_rejection(): + bob = BOB() + with pytest.raises(ValueError, match=r"m=0"): + bob.initialize_standalone(mf=MF, chif=CHIF, l=2, m=0) + + +# --------------------------------------------------------------------------- +# Test 6: w_r / tau overrides bypass Kerr lookup. +# --------------------------------------------------------------------------- + +def test_w_r_tau_overrides(): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF, l=2, m=2, w_r=0.5, tau=12.0) + assert bob.w_r == pytest.approx(0.5) + assert bob.tau == pytest.approx(12.0) + # Omega_QNM = w_r / |m|. + assert bob.Omega_QNM == pytest.approx(0.25) + + +# --------------------------------------------------------------------------- +# Test 7: construct_BOB short-circuit — no crash on self.data is None, +# and NR_based_on_BOB_ts stays None. Locks in the standalone branch in +# construct_BOB so a future refactor can't silently re-introduce the +# NR-alignment block unconditionally. +# --------------------------------------------------------------------------- + +def test_construct_BOB_short_circuit(): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF) + bob.what_should_BOB_create = "news" + # data is None in standalone mode — must not crash. + assert bob.data is None + t, y = bob.construct_BOB() + assert np.all(np.isfinite(y)) + # NR_based_on_BOB_ts is its default (None) — the alignment block didn't run. + assert bob.NR_based_on_BOB_ts is None + + +# --------------------------------------------------------------------------- +# Test 8: finite-t0 path works without NR. Only the optimize_* flags are +# rejected; the finite-t0 build itself is pure analytic math. +# --------------------------------------------------------------------------- + +def test_finite_t0_path(): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF) + # Order matters: set_initial_time requires what_should_BOB_create + # to be set first (existing API contract on the setter). + bob.what_should_BOB_create = "news" + omega0_before = bob.Omega_0 # captured at mode-set time + bob.set_initial_time = -50.0 # flips minf_t0 to False internally + # In standalone mode, set_initial_time skips the NR Omega_0 refit. + assert bob.Omega_0 == omega0_before + t, y = bob.construct_BOB() + assert np.all(np.isfinite(y)) + # The waveform was built on the finite-t0 path. + assert bob.minf_t0 is False + + +# --------------------------------------------------------------------------- +# Test 9: t0 + Omega_0 at init time. The user can pass both up front and +# then go straight to what_should_BOB_create + construct_BOB. +# --------------------------------------------------------------------------- + +def test_init_with_t0_and_Omega_0(): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF, t0=-50.0, Omega_0=0.123) + # Finite-t0 was selected automatically. + assert bob.minf_t0 is False + # t0 stored in absolute units (tp + relative); default tp=0 means t0=-50. + assert bob._wf_config.t0 == pytest.approx(-50.0) + # User Omega_0 stored and override flag set. + assert bob.Omega_0 == pytest.approx(0.123) + assert bob._runtime.Omega0_user_override is True + # User Omega_0 survives the mode-set step (the override gate in + # _apply_standalone_mode is the test target here). + bob.what_should_BOB_create = "news" + assert bob.Omega_0 == pytest.approx(0.123) + # End-to-end build works with no further setup beyond mode selection. + # Note: construct_BOB_finite_t0 resets Omega_0 to Omega_ISCO at the + # end as a documented cleanup hack (Stage 3.3 design decision), so we + # don't re-check Omega_0 after this call. + t, y = bob.construct_BOB() + assert np.all(np.isfinite(y)) + + +# --------------------------------------------------------------------------- +# Test 10: default Omega_0 applies the mode-appropriate fit at mode-set time. +# --------------------------------------------------------------------------- + +@pytest.mark.parametrize( + "mode, fit_fn", + [ + ("psi4", gen_utils.Omega_0_fit_psi4), + ("news", gen_utils.Omega_0_fit_news), + ("strain", gen_utils.Omega_0_fit_strain), + ], +) +def test_Omega_0_default_fit(mode, fit_fn): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF) + bob.what_should_BOB_create = mode + expected = fit_fn(bob.mf, bob.chif_with_sign) + # Bit-equality: omega_fn is a pure function called with identical inputs. + assert bob.Omega_0 == expected + + +def test_Omega_0_default_updates_on_mode_switch(): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF) + bob.what_should_BOB_create = "news" + after_news = bob.Omega_0 + bob.what_should_BOB_create = "psi4" + after_psi4 = bob.Omega_0 + # Mode switch refits Omega_0 for the new mode — they should differ. + assert after_news != after_psi4 + assert after_psi4 == gen_utils.Omega_0_fit_psi4(bob.mf, bob.chif_with_sign) + + +# --------------------------------------------------------------------------- +# Test 11: user-supplied Omega_0 is sticky across mode switches; the manual +# bob.Omega_0 setter still works. +# --------------------------------------------------------------------------- + +def test_Omega_0_user_override_sticky(): + bob = BOB() + bob.initialize_standalone(mf=MF, chif=CHIF, Omega_0=0.123) + bob.what_should_BOB_create = "news" + assert bob.Omega_0 == pytest.approx(0.123) + # Override flag is set. + assert bob._runtime.Omega0_user_override is True + + # Switch modes — override must persist. + bob.what_should_BOB_create = "psi4" + assert bob.Omega_0 == pytest.approx(0.123) + + # Manual setter still works for power users. + bob.Omega_0 = 0.5 + assert bob.Omega_0 == pytest.approx(0.5)