Normalize hybrid/sparse binary FB expansions and allow signed-H3 DDH fits

[vhaasteren]

Description

Two classes of published PTA timing models are mathematically well-defined but
currently rejected or incompletely represented by PINT (issue #2022 ).

1. Canonical orbital-frequency Taylor series

The binary orbital phase can be written as the number of elapsed orbits

N(dt) = FB0 dt
      + FB1 dt² / 2!
      + FB2 dt³ / 3!
      + ...
      + FBn dt^(n+1) / (n+1)! ,

where dt = t - Tref in seconds. The period representation supplies the first
two coefficients through

FB0 = 1 / (PB × 86400)
FB1 = -PBDOT × FB0².

Some published tempo2 par files use PB for the constant term while supplying
FB1, FB2, ... directly. PPTA DR2/DR3 J2241−5236, for example, provides
PB + FB1...FB17 but no FB0. Tempo2 evaluates that file by using
1/PB as the constant orbital frequency. PINT previously rejected it because
its FBX implementation required an explicit FB0.

This PR normalizes every such model to one canonical FBX representation:

• PB is converted exactly to the nominal FB0;
• an ordinary PBDOT is converted to FB1 when no explicit FB1 exists;
• explicit FB1 wins over redundant PBDOT, with a warning;
• missing interior Taylor coefficients are inserted as frozen zeros;
• PB and, when FB1 is present, PBDOT remain available as read-only
  derived views.

An explicitly valued ordinary PB together with an explicitly valued ordinary
FB0 remains an error and is rejected before normalization mutates the model.
For an existing pure-FBX model, this canonicalization is an API/output change:
PB and applicable PBDOT values become derived views,
as_parfile() gains commented # PB/# PBDOT lines, and an unset writable
PBDOT that the FBX orbit would ignore is removed.

Treating a missing interior coefficient as zero is important physically. For
example,

PB + FB2

means

FB0 = 1/(PB × 86400),   FB1 = 0,   FB2 = supplied value,

not “ignore FB2.” PINT's previous OrbitFBX implementation could silently
drop FB2+ after a missing index. This PR fixes that independent bug and adds
equation-level phase tests ensuring that every supplied coefficient actually
contributes.

Tempo2 ELL1 currently uses the presence of FB1 as a sentinel for the entire
higher-FB branch, so it ignores FB2+ when FB1 is omitted. PINT intentionally
uses the complete Taylor-series interpretation instead. The published
J2241−5236 model has a contiguous series and therefore remains exactly
tempo2-compatible.

The intentional timing deviations from current tempo2 are limited and
explicit:

• PB + FB2 and FB0 + FB2: PINT uses FB2 with an inserted FB1=0;
  tempo2 ELL1 currently ignores FB2;
• PB + PBDOT + FB2: PINT applies both the equivalent FB1 and FB2;
  tempo2 applies only the PBDOT-equivalent term;
• non-ELL1 use is a deliberate PINT generalization of its existing FBX phase
  machinery;
• DDGR+FBX is rejected by PINT, whereas current tempo2 DDGR silently uses
  PB and ignores FBn. The companion tempo2 PR should reject this
  unsupported combination too.

Conversely, FB0 + PBDOT has the same nominal phase in tempo2 ELL1 and in the
new PINT representation: PINT now converts it to FB1 instead of silently
ignoring it as current PINT does. A companion tempo2 issue/PR will replace the
FB1 sentinel with the same complete-series rule so the sparse cases converge
back to compatibility.

FBX is treated as the complete orbital-phase parameterization for every PINT
binary model that supports it, not as an additive correction to PB/PBDOT.
Accordingly, dynamically derived values such as DDGR's GR prediction for
PBDOT are never evaluated once and copied into an independent FB1; doing so
would become stale when masses change.

This PR explicitly rejects DDGR+FBX. A tempting normalization is to insert
FB1=0 and expose FB-derived PBDOT=0, but that fixes only the phase label:
DDGR derives post-Keplerian quantities directly from PB, and its analytic
derivatives are expressed with respect to PB. A physically correct FBX
implementation must keep those quantities synchronized through
PB=1/(86400 FB0) and add

d(PK)/d(FB0) = d(PK)/d(PB) × d(PB)/d(FB0),
d(PB)/d(FB0) = -1/(86400 FB0²).

Until that larger change is implemented and derivative-tested, accepting the
combination would produce an internally inconsistent timing model. A valued
XPBDOT with FBX is likewise rejected until a dynamic combination rule is
explicitly designed.

The nominal phase transformation is exact. Reported marginal uncertainties are
mapped with the first-order Jacobian:

sigma_FB0 = sigma_PB / (PB² × 86400)

sigma_FB1² =
    (FB0² sigma_PBDOT)²
  + (2 PBDOT FB0 sigma_FB0)².

The nonlinear change also creates FB0–FB1 covariance that a par file cannot
represent; this PR does not claim lossless covariance transfer.

2. DDH with a signed weak-signal H3 estimate

The DDH orthometric parameters map to the traditional Shapiro quantities as

SINI = 2 STIGMA / (1 + STIGMA²)
M2   = H3 / (Tsun STIGMA³).

For finite positive STIGMA, the Shapiro delay is

DeltaS = -2 H3 / STIGMA³ × log(B),

where B is the usual DD orbital geometry term. Wherever B>0, the delay is
linear and smooth in H3, including through H3=0; changing the sign of H3
does not change the logarithm's domain.

For a weak or absent Shapiro signal, an unconstrained amplitude estimate can
legitimately cross zero. MPTA DR2 J1825−0319 has such a best fit:

H3   < 0
STIG > 0.

The corresponding derived M2 is negative. That is not a physical companion
mass, but rejecting the model also rejects a finite, differentiable signed
amplitude fit that tempo2 can evaluate. It prevents faithful residual
comparison and can prevent a fitter from crossing the zero-signal point.

This PR therefore allows negative H3 in BinaryDDH and warns that the
negative derived M2 has no physical mass interpretation. The exception is
specific to DDH's smooth signed-H3 extension:

• free negative M2 still raises;
• generic derived SINI validation is unchanged;
• STIGMA must remain finite and positive;
• STIGMA=0 is rejected because the model contains inverse powers of it;
• other binary models receive no blanket funcParameter exemption.

Tempo2 evaluates negative H3 without PINT's derived-mass gate, so the nominal
DDH delay remains compatible. PINT intentionally remains stricter for
STIGMA<=0: negative values are outside the adopted physical domain and zero
is singular.

With the prototype normalization, the published PPTA DR3 J2241−5236 dataset
loads with finite residuals at approximately 0.821 µs RMS, and the MPTA DR2
J1825−0319 dataset loads as DDH with approximately 5.677 µs RMS, matching the
tempo2/libstempo timing evidence used to motivate the change.

Scope and tests

The implementation is confined to binary-parameter normalization, FBX gap
validation, DDH-native validation, documentation, and tests. Tests cover sparse
series, phase equivalence, uncertainties, fit flags, serialization,
idempotency, TCB conversion, general binary models, DDGR+FBX rejection before
mutation, negative-H3 continuity and fitting, and unchanged free-M2/SINI
validation.

Reported by Wang-Wei Yu @Wang-weiYu during IPTA DR3 / 5-PTA combination testing with
the updated MetaPulsar version.