Vectorized non-GPU time marginalization offsets by ONE detector's max while exponentiating the accumulated network sum (platform-conditional; float64 platforms get #232's 745-nat budget)

[oshaughnessy-junior]

Found while sweeping the class behind #232. This is a different defect — same family (a mis-taken offset before exp), but per-row rather than batch, and it is not currently producing wrong numbers on the production platform. Filing it rather than patching it, because its impact is platform-conditional and it deserves its own reachability call.

The site

MonteCarloMarginalizeCode/Code/RIFT/likelihood/factored_likelihood.py, in DiscreteFactoredLogLikelihoodViaArrayVector (def ~1992 on rift_O4d; the reachable --vectorized path without --gpu), inside the for det in detectors loop:

            lnL = term1+term2
            lnL_array[indx_ex] += lnL      # accumulates over detectors
            maxlnL = np.max(lnL)           # <-- only THIS detector's contribution
            lnLmargOut[indx_ex] = maxlnL + np.log(my_simps(np.exp(lnL_array[indx_ex] - maxlnL), dx=deltaT))

The offset maxlnL is the max of the current detector's own term, while the array being exponentiated, lnL_array[indx_ex], is the running sum over all detectors so far. The statement also sits inside the detector loop, so it is recomputed per detector and only the last iteration survives — correct in exact arithmetic (the last pass holds the full sum), but the surviving offset is the last detector's max applied to the whole network's sum.

The expression is algebraically offset-invariant, so this is a floating-point range defect, not an arithmetic error. For an N-detector network the argument to exp sits roughly (N−1) offsets too high.

Why it is not biting today, and where it would

lnL_array is RiftFloat. On x86_64 Linux that is 80-bit longdouble:

platform	RiftFloat	exp underflow budget
x86_64 Linux (CIT, OSG)	longdouble	11 355 nats
macOS arm64, Windows MSVC, non-x86 Linux	float64	745 nats

That second row is not speculation — RIFT/precision.py documents exactly this fallback and even exports RIFT_FLOAT_HIGH_PRECISION for code that needs to know. So on any platform where RIFT_FLOAT_HIGH_PRECISION is False, this site carries the same 745-nat budget as #232, on a path that #232's fix does not touch.

At production amplitudes on x86_64 the margin is comfortable and I have no measurement of a wrong number from it. I am explicitly not claiming one.

What it does not explain

While measuring #232 I saw this path return -inf for ~0.7 % of draws (71 per 10 000 at cycle 1, rising to ~13 % as the sampler contracted). That is not this defect: the quadratic self-term is -(rho_h^2/2)(d_ref/d)^2, so a draw at the small-distance end of the prior sits of order -3e7 nats and underflows in any precision, 80-bit included. Those -inf are physical, and the AV sampler handles them. Recording it so the two are not conflated later.

Suggested fix, if wanted

Offset by the accumulated row, and hoist the reduction out of the detector loop so it runs once on the completed sum:

    # after the detector loop
    maxlnL = np.max(lnL_array, axis=-1, keepdims=True)
    lnLmargOut = maxlnL[...,0] + np.log(my_simps(np.exp(lnL_array - maxlnL), dx=deltaT, axis=-1))

That is the idiom already used at factored_likelihood_with_rotation.py:948, factored_likelihood_freqresponse.py:457, jax_ile/core.py, and now (via #187/#234) at the #232 site. It also removes N−1 redundant Simpson integrations per call.

Caveat for whoever takes it: it changes float64 rounding on a reachable production path, so it needs the same treatment #234 got — a test watched failing first, and a measured low-amplitude agreement number rather than a bit-identity claim. It should probably be reproduced on a RIFT_FLOAT_HIGH_PRECISION == False platform, since that is where it actually matters.

Related

• Vectorized time-marginalized likelihood offsets by the BATCH max, not the per-sample max: lnL underflows to -inf above rho ~ 40 and mcsamplerAV collapses (device-independent; affects a delivered O4b result) #232 — the batch-max defect, fixed by factored_likelihood: offset the vectorized time integral per row, not per batch (fixes -inf / mcsamplerAV collapse above rho ~ 40; issue #232) oshaughn/research-projects-RIT#187 (rift_O4c) and factored_likelihood: offset the vectorized time marginalization per extrinsic sample, not per batch (#232) #234 (rift_O4d).
• Also swept and found correct (genuine 1-D per-sample reductions): factored_likelihood.py ~828, ~1988, ~2074-equivalents on both branches; with_rotation, freqresponse, jax_ile, LISA.
• DiscreteFactoredLogLikelihoodViaArrayVectorNoLoopOrig carries the Vectorized time-marginalized likelihood offsets by the BATCH max, not the per-sample max: lnL underflows to -inf above rho ~ 40 and mcsamplerAV collapses (device-independent; affects a delivered O4b result) #232 construction verbatim but has zero production callers on either branch — the repo's own instrumented probe records {'...NoLoopOrig': 0} in RIFT/integrators/VALIDATION_rvs_weight_migration.md. Left unpatched deliberately; the one-line fix is the same as D1: AlternateIteration's extrinsic stage read the previous iteration's grid #187's if it is ever revived.