ported/financial_time_value_family

fn factor(self, rate: f64) -> f64 {
    match self {
        Self::EndOfPeriod => 1.0,
        Self::BeginningOfPeriod => 1.0 + rate,
    }
}

fn validate_finite(values: &[f64]) -> Result<(), FinancialError> {
    if values.iter().all(|value| value.is_finite()) {
        Ok(())
    } else {
        Err(FinancialError::Value)
    }
}

fn growth(periodic_rate: f64, periods: f64) -> Result<f64, FinancialError> {
    validate_finite(&[periodic_rate, periods])?;
    let base = 1.0 + periodic_rate;
    if base <= 0.0 {
        return Err(FinancialError::Num);
    }
    let factor = match power_kernel(base, periods) {
        Ok(factor) => factor,
        Err(WorksheetErrorCode::Div0) => return Err(FinancialError::Div0),
        Err(WorksheetErrorCode::Num) => return Err(FinancialError::Num),
        Err(_) => return Err(FinancialError::Value),
    };
    if factor.is_finite() {
        Ok(factor)
    } else {
        Err(FinancialError::Num)
    }
}

fn annuity_term(
    periodic_rate: f64,
    periods: f64,
    timing: PaymentTiming,
) -> Result<f64, FinancialError> {
    if periodic_rate.abs() < EPSILON {
        return Ok(periods);
    }
    let factor = growth(periodic_rate, periods)?;
    let term = timing.factor(periodic_rate) * (factor - 1.0) / periodic_rate;
    if term.is_finite() {
        Ok(term)
    } else {
        Err(FinancialError::Num)
    }
}

fn growth_from_base(base: f64, periods: f64) -> Result<f64, FinancialError> {
    let factor = match power_kernel(base, periods) {
        Ok(factor) => factor,
        Err(WorksheetErrorCode::Div0) => return Err(FinancialError::Div0),
        Err(WorksheetErrorCode::Num) => return Err(FinancialError::Num),
        Err(_) => return Err(FinancialError::Value),
    };
    if factor.is_finite() {
        Ok(factor)
    } else {
        Err(FinancialError::Num)
    }
}

fn balance_equation(
    periodic_rate: f64,
    periods: f64,
    payment_value: f64,
    present_value: f64,
    future_value: f64,
    timing: PaymentTiming,
) -> Result<f64, FinancialError> {
    validate_finite(&[
        periodic_rate,
        periods,
        payment_value,
        present_value,
        future_value,
    ])?;
    if periodic_rate.abs() < EPSILON {
        return Ok(future_value + present_value + payment_value * periods);
    }
    let factor = growth(periodic_rate, periods)?;
    let term = annuity_term(periodic_rate, periods, timing)?;
    let value = future_value + present_value * factor + payment_value * term;
    if value.is_finite() {
        Ok(value)
    } else {
        Err(FinancialError::Num)
    }
}

/// PMT on the correctly-rounded numeric core (W108 Bead C).
///
/// PV-side stable discount form (the live-Excel search winner). Instead of the
/// FV `powi` factor + `(F-1)/r` term (catastrophically wrong on the PMT lane —
/// up to ~5.5e8 ULP), this reproduces Excel's own `exp`/`log1p`/`expm1` chain:
///
/// ```text
/// negL  = -(n * log1p(r))        // -n·ln(1+r)
/// invF  = exp(negL)              // (1+r)^-n   (no cancellation)
/// denom = -expm1(negL)           // 1 - (1+r)^-n  (stable; avoids 1 - tiny)
/// tf    = 1 + r·type
/// pmt   = -(pv + fv·invF) · r / (tf · denom)
/// ```
///
/// The transcendentals are the faithful (~0.5-ULP) deterministic core, matching
/// Excel's bespoke routine on ~99% of inputs. The residual is the tracked W108
/// floor (≤2-3 ULP open discrepancy), NOT an acceptance.
pub fn pmt(
    periodic_rate: f64,
    periods: f64,
    present_value: f64,
    future_value: f64,
    timing: PaymentTiming,
) -> Result<f64, FinancialError> {
    validate_finite(&[periodic_rate, periods, present_value, future_value])?;
    if periods == 0.0 {
        return Err(FinancialError::Num);
    }
    // W108 Bead C: the zero-rate special case is EXACT `r == 0.0`, matching Excel
    // and the live-search model — a merely-tiny rate (e.g. 1e-13) is NOT collapsed
    // to the linear branch (that was a ~3813-ULP miss); it flows through the
    // stable discount form (1-ULP on PMT-0025).
    if periodic_rate == 0.0 {
        return Ok(-(future_value + present_value) / periods);
    }
    // Growth-domain guard: (1+r) <= 0 has no real (1+r)^-n — #NUM!, matching the
    // prior `growth` guard and live Excel.
    if 1.0 + periodic_rate <= 0.0 {
        return Err(FinancialError::Num);
    }
    // R31 adjacent-rate calls: beginning payments reject the pole in the
    // stored-PC64 reciprocal coefficient, even when balances cancel. End
    // payments at the same rate retain their ordinary discount domain.
    if matches!(timing, PaymentTiming::BeginningOfPeriod)
        && excel_x87_add(1.0, excel_x87_recip(periodic_rate)) == 0.0
    {
        return Err(FinancialError::Num);
    }
    // W109 G6-01 (2026-07-23, multi-agent workflow): the transcendental substrate is
    // the x87 chain (Excel's PMT `exp` is the x87 exponential `excel_exp`, matched 234/234; `expm1` is the internal
    // x87 double-rounded Kahan), and the combine is the live-Excel-pinned **quotient-first,
    // ·rate LAST** graph for end payments. R30 separately identifies the positive-rate
    // discount/type-one coefficient and its stores below; other routes keep the timing division.
    //
    // dc-j9k.12: the logarithm is not a correctly rounded log1p but the atanh series of
    // `excel_log1p_atanh`. Holding pv varied at fixed (rate, nper) pins Excel's binary64 em in
    // 4,720 recorded groups (fv 0, type 0, rate > 0). For each of the 1,184 rates among them
    // one binary64 L reproduces em = expm1(-(nper * L)) through the x87 Kahan expm1 below at
    // every nper; the series gives that L for 1,181 rates, a correctly rounded log1p for 667.
    // The series is used where its twelve terms have converged (-1/3 <= rate <= 1/2, so
    // |s| <= 1/5); outside that range Excel's logarithm is not yet identified and the
    // correctly rounded one stays.
    let log_growth = if (-1.0 / 3.0..=0.5).contains(&periodic_rate) {
        excel_log1p_atanh(periodic_rate)
    } else {
        excel_log1p(periodic_rate)
    };
    let neg_log = -(periods * log_growth);
    let em = crate::ported::excel_numeric::excel_expm1_spilled(neg_log); // (1+r)^-n − 1  (x87 Kahan; negative for r>0)
    // Preserve the growth-reciprocal overflow refusal before the balance
    // fold can turn pv/Inf + fv into a finite payment. R31 paired parts keep
    // first-period beginning IPMT zero ahead of this shared PMT check.
    if em == 0.0 || !em.is_finite() {
        return Err(FinancialError::Num);
    }
    let tf = timing.factor(periodic_rate); // 1 + r·type
    // Equal opposite balances are interest only: pv + fv*v = fv*(v - 1).
    // Cancelling that common factor before rounding avoids subtracting nearly
    // equal balances. Negative rates use the common stored fold below: its
    // quotient is exactly fv here, retaining their measured final stores.
    if periodic_rate > 0.0 && present_value != 0.0 && future_value == -present_value {
        let product = future_value * periodic_rate;
        let result = if product.is_finite() {
            product / tf
        } else {
            let scaled = (future_value / tf) * periodic_rate;
            if scaled.is_finite() {
                scaled
            } else {
                // With beginning payments and a positive rate, rate/tf <= 1.
                // This order cannot overflow a finite balance. The previous
                // order can round MAX/tf upward and overflow its final product.
                future_value * (periodic_rate / tf)
            }
        };
        return if result.is_finite() { Ok(result) } else { Err(FinancialError::Num) };
    }
    // Both arms store the quotient from PC64 (dc-j9k.12). With pv varied at fixed (rate, nper),
    // fv 0 and type 0, q = RN53(RN64(pv / em)) followed by RN53(RN64(q * rate)) admits one
    // binary64 em in every one of 4,802 recorded groups; a binary64 quotient leaves 59 groups
    // with no em at all.
    let quotient = if neg_log < 0.0 || periodic_rate < 0.0 {
        // dc-j9k.12, R29/R31: the discount and negative-rate arms combine balances as
        // (pv + fv)/em + fv, with both additions and the division stored
        // from PC64. Public type-zero calls constrain em independently:
        // this fold matches all192 paired future-value calls for each rate sign.
        // Separate negative-rate ties distinguish each addition (84 and114 readings).
        // Forming pv + fv*(1+em) loses those recorded roundings.
        let balance = excel_x87_add(present_value, future_value);
        excel_x87_add(excel_x87_div(balance, em), future_value)
    } else {
        let inv_factor = 1.0 + em;
        let numerator = present_value + future_value * inv_factor;
        excel_x87_div(numerator, em)
    };
    let result = if (neg_log < 0.0 || periodic_rate < 0.0)
        && matches!(timing, PaymentTiming::BeginningOfPeriod)
    {
        // R30/R31: paired type-zero calls pin em independently of timing. The
        // discount and negative-rate arms multiply by 1/(1+1/r), storing each
        // operation from PC64. Negative-rate calls also pin this coefficient
        // under either final-word model before 44 fresh product ties distinguish
        // the PC64 store. Dividing by (1+r) loses those recorded roundings.
        let reciprocal = excel_x87_recip(periodic_rate);
        let factor = if reciprocal.is_finite() {
            excel_x87_recip(excel_x87_add(1.0, reciprocal))
        } else {
            // A rate whose reciprocal overflows is subnormal. Its exact
            // r/(1+r) rounds to r: r^2 is far below half a subnormal unit.
            periodic_rate
        };
        excel_x87_mul(quotient, factor)
    } else {
        let quotient = excel_x87_div(quotient, tf);
        // R31: 32 independently measured negative-rate product ties require
        // the PC64-to-binary64 store after the already stored quotient.
        // Their denominators are pinned by separate public type-zero calls.
        if neg_log < 0.0 || periodic_rate < 0.0 {
            excel_x87_mul(quotient, periodic_rate)
        } else {
            quotient * periodic_rate
        }
    };
    if result.is_finite() {
        Ok(result)
    } else {
        Err(FinancialError::Num)
    }
}

fn values_bracket_root(lower_value: f64, upper_value: f64) -> bool {
    (lower_value < 0.0 && upper_value > 0.0) || (lower_value > 0.0 && upper_value < 0.0)
}

fn bisect_rate_bracket(
    mut lower: f64,
    mut lower_value: f64,
    mut upper: f64,
    upper_value: f64,
    periods: f64,
    payment_value: f64,
    present_value: f64,
    future_value: f64,
    timing: PaymentTiming,
) -> Result<f64, FinancialError> {
    if lower_value.abs() <= RATE_TOLERANCE {
        return Ok(lower);
    }
    if upper_value.abs() <= RATE_TOLERANCE {
        return Ok(upper);
    }
    if !values_bracket_root(lower_value, upper_value) {
        return Err(FinancialError::NoConvergence);
    }

    for _ in 0..RATE_BRACKET_MAX_ITERATIONS {
        let midpoint = (lower + upper) / 2.0;
        let midpoint_value = balance_equation(
            midpoint,
            periods,
            payment_value,
            present_value,
            future_value,
            timing,
        )?;
        if midpoint_value.abs() <= RATE_TOLERANCE {
            return Ok(midpoint);
        }

        if values_bracket_root(lower_value, midpoint_value) {
            upper = midpoint;
        } else {
            lower = midpoint;
            lower_value = midpoint_value;
        }
    }

    let midpoint = (lower + upper) / 2.0;
    let midpoint_value = balance_equation(
        midpoint,
        periods,
        payment_value,
        present_value,
        future_value,
        timing,
    )?;
    if midpoint_value.abs() <= RATE_TOLERANCE || (upper - lower).abs() <= EPSILON {
        Ok(midpoint)
    } else {
        Err(FinancialError::NoConvergence)
    }
}

fn search_rate_bracket_from_zero(
    zero_value: f64,
    initial_rate: f64,
    periods: f64,
    payment_value: f64,
    present_value: f64,
    future_value: f64,
    timing: PaymentTiming,
) -> Result<Option<f64>, FinancialError> {
    let mut positive_rate = initial_rate.abs().max(0.1);
    for _ in 0..RATE_BRACKET_MAX_ITERATIONS {
        let positive_value = balance_equation(
            positive_rate,
            periods,
            payment_value,
            present_value,
            future_value,
            timing,
        )?;
        if positive_value.abs() <= RATE_TOLERANCE {
            return Ok(Some(positive_rate));
        }
        if values_bracket_root(zero_value, positive_value) {
            return bisect_rate_bracket(
                0.0,
                zero_value,
                positive_rate,
                positive_value,
                periods,
                payment_value,
                present_value,
                future_value,
                timing,
            )
            .map(Some);
        }

        positive_rate *= 2.0;
        if !positive_rate.is_finite() {
            break;
        }
    }

    let mut negative_rate = if initial_rate < 0.0 {
        initial_rate
    } else {
        -0.1
    };
    if negative_rate <= MIN_VALID_RATE {
        negative_rate = -0.1;
    }
    for _ in 0..RATE_BRACKET_MAX_ITERATIONS {
        let negative_value = balance_equation(
            negative_rate,
            periods,
            payment_value,
            present_value,
            future_value,
            timing,
        )?;
        if negative_value.abs() <= RATE_TOLERANCE {
            return Ok(Some(negative_rate));
        }
        if values_bracket_root(negative_value, zero_value) {
            return bisect_rate_bracket(
                negative_rate,
                negative_value,
                0.0,
                zero_value,
                periods,
                payment_value,
                present_value,
                future_value,
                timing,
            )
            .map(Some);
        }

        let margin = 1.0 + negative_rate;
        if margin <= EPSILON {
            break;
        }
        negative_rate = -1.0 + margin / 2.0;
        if negative_rate <= MIN_VALID_RATE {
            negative_rate = MIN_VALID_RATE + EPSILON;
        }
    }

    Ok(None)
}

pub fn rate(
    periods: f64,
    payment_value: f64,
    present_value: f64,
    future_value: f64,
    timing: PaymentTiming,
    guess: Option<f64>,
) -> Result<f64, FinancialError> {
    let root = rate_root(periods, payment_value, present_value, future_value, timing, guess)?;
    Ok(polish_root(root, |r| {
        balance_equation(r, periods, payment_value, present_value, future_value, timing)
    }))
}

/// A root `seed` of `f`, polished by Newton steps (a central-difference slope) to full binary64
/// convergence (dc-ifc): Excel's RATE converges fully where the 1e-7 tolerance of the search
/// stops early (RATE(48,-200,8000): Excel 7e-16 from the exact root, the search's root 2e-13).
/// The steps stop when they no longer contract, keeping whichever adjacent value has the smaller
/// residual, or when `f` fails, keeping the last value.
pub(crate) fn polish_root(seed: f64, f: impl Fn(f64) -> Result<f64, FinancialError>) -> f64 {
    let mut rate = seed;
    let mut last_step = f64::INFINITY;
    for _ in 0..60 {
        let Ok(value) = f(rate) else { break };
        if value == 0.0 {
            break;
        }
        let h = rate.abs().max(1e-4) * 1e-6;
        let (Ok(up), Ok(down)) = (f(rate + h), f(rate - h)) else { break };
        let slope = (up - down) / (2.0 * h);
        if !slope.is_finite() || slope == 0.0 {
            break;
        }
        let next = rate - value / slope;
        if !next.is_finite() || next.to_bits() == rate.to_bits() || next <= MIN_VALID_RATE {
            break;
        }
        let step = (next - rate).abs();
        if step >= last_step {
            if f(next).is_ok_and(|v| v.abs() <= value.abs()) {
                rate = next;
            }
            break;
        }
        last_step = step;
        rate = next;
    }
    rate
}

/// RATE's root search: secant steps from the guess, bisection once a bracket is found, to a
/// residual of at most [`RATE_TOLERANCE`].
fn rate_root(
    periods: f64,
    payment_value: f64,
    present_value: f64,
    future_value: f64,
    timing: PaymentTiming,
    guess: Option<f64>,
) -> Result<f64, FinancialError> {
    validate_finite(&[periods, payment_value, present_value, future_value])?;
    let mut prev_rate = 0.0;
    let mut prev_value = balance_equation(
        prev_rate,
        periods,
        payment_value,
        present_value,
        future_value,
        timing,
    )?;
    let mut current_rate = guess.unwrap_or(0.1);
    if current_rate <= MIN_VALID_RATE {
        current_rate = MIN_VALID_RATE;
    }
    let mut current_value = balance_equation(
        current_rate,
        periods,
        payment_value,
        present_value,
        future_value,
        timing,
    )?;

    if current_value.abs() <= RATE_TOLERANCE {
        return Ok(current_rate);
    }
    if values_bracket_root(prev_value, current_value) {
        return bisect_rate_bracket(
            prev_rate,
            prev_value,
            current_rate,
            current_value,
            periods,
            payment_value,
            present_value,
            future_value,
            timing,
        );
    }

    for _ in 0..RATE_MAX_ITERATIONS {
        if current_value.abs() <= RATE_TOLERANCE {
            return Ok(current_rate);
        }

        let denominator = current_value - prev_value;
        let next_rate = if denominator.abs() > EPSILON {
            current_rate - current_value * (current_rate - prev_rate) / denominator
        } else {
            let h = current_rate.abs().max(1.0) * 1e-6;
            let lower = (current_rate - h).max(-0.999_999_999);
            let upper = current_rate + h;
            let lower_value = balance_equation(
                lower,
                periods,
                payment_value,
                present_value,
                future_value,
                timing,
            )?;
            let upper_value = balance_equation(
                upper,
                periods,
                payment_value,
                present_value,
                future_value,
                timing,
            )?;
            let derivative = (upper_value - lower_value) / (upper - lower);
            if derivative.abs() <= EPSILON {
                return Err(FinancialError::NoConvergence);
            }
            current_rate - current_value / derivative
        };

        if !next_rate.is_finite() || next_rate <= MIN_VALID_RATE {
            break;
        }

        prev_rate = current_rate;
        prev_value = current_value;
        current_rate = next_rate;
        current_value = balance_equation(
            current_rate,
            periods,
            payment_value,
            present_value,
            future_value,
            timing,
        )?;

        if values_bracket_root(prev_value, current_value) {
            return bisect_rate_bracket(
                prev_rate,
                prev_value,
                current_rate,
                current_value,
                periods,
                payment_value,
                present_value,
                future_value,
                timing,
            );
        }
    }

    if current_value.abs() <= RATE_TOLERANCE {
        Ok(current_rate)
    } else if let Some(root) = search_rate_bracket_from_zero(
        balance_equation(
            0.0,
            periods,
            payment_value,
            present_value,
            future_value,
            timing,
        )?,
        guess.unwrap_or(0.1),
        periods,
        payment_value,
        present_value,
        future_value,
        timing,
    )? {
        Ok(root)
    } else {
        Err(FinancialError::NoConvergence)
    }
}

/// What Excel's PMT, IPMT and PPMT share for payments at the end of a period (dc-j9k.12): the
/// series logarithm `L` of 1 + rate and PMT's stored quotient
/// `q = RN53(RN64((pv + fv * v^n) / em_n))`, `em_n = expm1(-(n * L))`, `v^n = 1 + em_n`.
/// `None` where PMT itself takes another route: a rate outside the series' range, equal and
/// opposite balances, or an em of zero.
fn end_payment_quotient(
    periodic_rate: f64,
    periods: f64,
    present_value: f64,
    future_value: f64,
) -> Option<(f64, f64)> {
    if periodic_rate == 0.0 || !(-1.0 / 3.0..=0.5).contains(&periodic_rate) {
        return None;
    }
    if present_value != 0.0 && future_value == -present_value {
        return None;
    }
    let log_growth = excel_log1p_atanh(periodic_rate);
    let em = crate::ported::excel_numeric::excel_expm1_spilled(-(periods * log_growth));
    if em == 0.0 || !em.is_finite() {
        return None;
    }
    let quotient = excel_x87_div(present_value + future_value * (1.0 + em), em);
    quotient.is_finite().then_some((log_growth, quotient))
}

/// Period `period_index` of an end-of-period annuity on Excel's own pieces (dc-j9k.12), with
/// `m = periods - period_index + 1` periods left, `t = -(m * L)`, `em_m = expm1(t)` (the x87
/// Kahan expm1) and `u = exp(t)` (the x87 exponential that expm1 starts from):
/// interest `-RN(RN(RN(q * em_m) - RN(fv * u)) * rate)` and principal
/// `RN(RN(RN(q - fv) * u) * rate)`, every RN the x87 store RN53(RN64(.)).
fn end_payment_parts(
    periodic_rate: f64,
    period_index: f64,
    periods: f64,
    future_value: f64,
    log_growth: f64,
    quotient: f64,
) -> (f64, f64) {
    let t = -((periods - period_index + 1.0) * log_growth);
    let em = crate::ported::excel_numeric::excel_expm1_spilled(t);
    let u = crate::ported::excel_numeric::excel_exp(t);
    let interest = -excel_x87_mul(
        excel_x87_sub(excel_x87_mul(quotient, em), excel_x87_mul(future_value, u)),
        periodic_rate,
    );
    let principal = excel_x87_mul(excel_x87_mul(excel_x87_sub(quotient, future_value), u), periodic_rate);
    (interest, principal)
}

/// IPMT on the correctly-rounded numeric core (W108 Bead C). Interest portion of
/// payment `period_index`: the remaining balance entering the period times `r`,
/// reconstructed from the core PMT and the `fv_exp` balance recurrence.
pub fn ipmt(
    periodic_rate: f64,
    period_index: f64,
    periods: f64,
    present_value: f64,
    future_value: f64,
    timing: PaymentTiming,
) -> Result<f64, FinancialError> {
    validate_finite(&[
        periodic_rate,
        period_index,
        periods,
        present_value,
        future_value,
    ])?;
    if period_index < 1.0 || period_index > periods {
        return Err(FinancialError::Num);
    }
    // R31 public IPMT calls return zero here even at rate -1 and below.
    // This precedes PMT's rate and reciprocal-pole checks, but follows the
    // finite-argument and period-range checks.
    if period_index == 1.0 && matches!(timing, PaymentTiming::BeginningOfPeriod) {
        return Ok(0.0);
    }
    // W109 G6-01 (2026-07-23): EXACT zero-rate special case (mirrors PMT's W108 Bead C),
    // NOT an `abs() < EPSILON` band. Live-Excel: IPMT(0,per,..)=+0.0, but a merely-tiny rate
    // (e.g. ±1e-13) flows through the MAIN path (IPMT(1e-13,1,360,2e5)=−2e-8, IPMT(−1e-13,..)=
    // +2e-8) — the old band wrongly collapsed tiny/negative rates to 0.0 and dropped the sign.
    if periodic_rate == 0.0 {
        return Ok(0.0);
    }
    if 1.0 + periodic_rate <= 0.0 {
        return Err(FinancialError::Num);
    }
    let periodic_payment = pmt(periodic_rate, periods, present_value, future_value, timing)?;
    // dc-j9k.12: payments at the end of a period take Excel's own pieces. On the catalogue
    // oracle and bank rows with fv 0 (95 of 95) the interest is -RN(RN(q * em_m) * rate).
    if matches!(timing, PaymentTiming::EndOfPeriod)
        && let Some((log_growth, quotient)) =
            end_payment_quotient(periodic_rate, periods, present_value, future_value)
    {
        let (interest, _) = end_payment_parts(
            periodic_rate, period_index, periods, future_value, log_growth, quotient,
        );
        return if interest.is_finite() { Ok(interest) } else { Err(FinancialError::Num) };
    }
    let tf = timing.factor(periodic_rate); // 1 + r·type
    let type_is_beginning = matches!(timing, PaymentTiming::BeginningOfPeriod);
    // Remaining balance entering `period_index` (the `rip` term); interest = rip·r. Excel
    // computes it backward from the end of the annuity (dc-j9k): the balance B after
    // `period_index - 1` periods is the one that `m = periods - period_index + 1` more periods
    // take to -fv, B = -fv·v^m + pmt·(1 + r·type)·(v^m - 1)/r with v = 1/(1 + r), so late
    // periods of a long, high-rate annuity do not lose the balance to cancellation. Against
    // the catalogue oracle's IPMT and PPMT readings, the forward recurrence missed by up to
    // 0.6 in absolute terms (PPMT(0.1713,191,264,-8172.57) was -0.588 where Excel gives
    // 0.0116); the backward form is within a few ulp of Excel on every reading.
    let rip = if period_index <= 1.0 {
        if type_is_beginning {
            0.0
        } else {
            -present_value
        }
    } else {
        let remaining = periods - period_index + 1.0;
        let log_discount = -(remaining * excel_log1p(periodic_rate));
        let discount = exp_portable(log_discount); // v^m
        let discount_m1 = excel_expm1(log_discount); // v^m - 1
        let balance =
            -future_value * discount + periodic_payment * tf * (discount_m1 / periodic_rate);
        if type_is_beginning {
            -balance / (1.0 + periodic_rate)
        } else {
            -balance
        }
    };
    let result = rip * periodic_rate;
    if result.is_finite() {
        Ok(result)
    } else {
        Err(FinancialError::Num)
    }
}

pub fn ppmt(
    periodic_rate: f64,
    period_index: f64,
    periods: f64,
    present_value: f64,
    future_value: f64,
    timing: PaymentTiming,
) -> Result<f64, FinancialError> {
    let payment = pmt(periodic_rate, periods, present_value, future_value, timing)?;
    // IPMT's checks (period range, rate) and its zero-rate case apply to PPMT too.
    let interest = ipmt(
        periodic_rate,
        period_index,
        periods,
        present_value,
        future_value,
        timing,
    )?;
    if periodic_rate == 0.0 || (period_index <= 1.0 && matches!(timing, PaymentTiming::BeginningOfPeriod)) {
        return Ok(payment - interest);
    }
    // dc-j9k.12: payments at the end of a period take Excel's own pieces: the principal is
    // RN(RN(RN(q - fv) * u) * rate) with u = exp(-(m * L)), not a payment less interest.
    // Catalogue oracle, bank and fresh Excel rows: 138 of 138 with fv 0, 25 of 27 with fv.
    if matches!(timing, PaymentTiming::EndOfPeriod)
        && let Some((log_growth, quotient)) =
            end_payment_quotient(periodic_rate, periods, present_value, future_value)
    {
        let (_, principal) = end_payment_parts(
            periodic_rate, period_index, periods, future_value, log_growth, quotient,
        );
        return if principal.is_finite() { Ok(principal) } else { Err(FinancialError::Num) };
    }
    // Excel does not subtract the interest from the payment, which cancels when the interest
    // is nearly the whole payment (dc-j9k: PPMT(0.1809,1,375,57755.66,635.76) is -8.79e-24 in
    // Excel, and payment - interest gave noise of order 1e-12). The principal is the payment
    // net of the future value's interest, discounted over the m = periods - period_index + 1
    // periods left: v^m * (pmt - fv*r), and v^m * (pmt - fv*r/(1 + r)) for payments at the
    // start of a period. The catalogue oracle's PPMT readings hold it.
    let remaining = periods - period_index + 1.0;
    let discount = exp_portable(-(remaining * excel_log1p(periodic_rate)));
    let carried = if matches!(timing, PaymentTiming::BeginningOfPeriod) {
        future_value * periodic_rate / (1.0 + periodic_rate)
    } else {
        future_value * periodic_rate
    };
    let result = discount * (payment - carried);
    if result.is_finite() {
        Ok(result)
    } else {
        Err(FinancialError::Num)
    }
}