"""PREREG_estimated_transform_2026-09-22 (control K3): Chandrupatla's method on
g = sign(f)|f|^(1/m_hat), with m_hat estimated during the run.

The solver is simroot.baselines.chandrupatla line for line (it matches scipy's
elementwise.find_root on all 199 problems); it keeps the TRUE values of f at its three points, so
signs, the best point and the stopping rule are those of f, and only the acceptance band's
(xi, phi) and the inverse quadratic interpolation read G(v) = sign(v)|v|^(1/m_hat), recomputed
from the stored f whenever m_hat changes - no extra evaluation.  m_hat starts at 1 and is
replaced by an order read from the last three evaluations only when it lies on a branch of
section 6.1 (m >= 2 or 0 < m <= 2/3); there is no reset.  With estimator="none" the method is
the original, bit for bit (tested).
"""
from __future__ import annotations

import math

from .baselines import Counter, _Budget
from .pure import _analogy_order, _config3, _lagouanelle_raw

# "lag2": POSTHOC_lag_branch_2026-09-22 - Lagouanelle on the branch m >= 2 only (post hoc)
ESTIMATORS = ("none", "lag", "analogy", "lag2")


def _G(v, inv):
    if inv == 1.0:
        return v
    if v == 0.0:
        return 0.0
    try:
        return math.copysign(abs(v) ** inv, v)
    except OverflowError:
        # 1/m_hat > 1 on the fractional branch: a huge |f| overflows; the interpolation that
        # reads it is then not finite and the step halves, as the method does for any such t
        return math.copysign(math.inf, v)


def chandrupatla_transformed(f, a, b, xtol, maxfev, estimator="none", stats=None):
    if estimator not in ESTIMATORS:
        raise ValueError(f"unknown estimator {estimator!r}")
    g = Counter(f, maxfev)
    x1, x2 = b, a
    f1 = f2 = 0.0
    m_hat = 1.0
    hist = []             # evaluations in order: (x, f)
    configs = []          # analogy memory (estimator "analogy")

    def observe(x, fx):
        nonlocal m_hat
        hist.append((x, fx))
        if estimator == "none" or len(hist) < 3:
            return
        (xp2, fp2), (xp, fp), (xc, fc) = hist[-3], hist[-2], hist[-1]
        if estimator in ("lag", "lag2"):
            m = _lagouanelle_raw(xc, fc, xp, fp, xp2, fp2)
            if estimator == "lag2" and m is not None and m < 2.0:
                m = None            # the fractional branch is not used (post hoc variant)
        else:
            cur = _config3(xc, fc, xp, fp, xp2, fp2)
            m = _analogy_order(cur, configs, translation=True)
            if cur is not None:
                configs.append(cur)
        if m is not None and (m >= 2.0 or 0.0 < m <= 2.0 / 3.0):
            if stats is not None and m != m_hat:
                stats["updates"] = stats.get("updates", 0) + 1
            m_hat = m
            if stats is not None:
                stats["m_hat"] = m_hat      # diagnostic only: the order last held

    try:
        f1 = g(x1)
        observe(x1, f1)
        f2 = g(x2)
        observe(x2, f2)
        if f1 == 0.0:
            return x1, g.n, True
        if f2 == 0.0:
            return x2, g.n, True
        x3, f3 = x2, f2
        t = 0.5
        while True:
            x = x1 + t * (x2 - x1)
            fx = g(x)
            observe(x, fx)
            if (fx > 0.0) == (f1 > 0.0):
                x3, f3 = x1, f1
                x1, f1 = x, fx
            else:
                x3, f3 = x2, f2
                x2, f2 = x1, f1
                x1, f1 = x, fx
            xm, fm = (x1, f1) if abs(f1) < abs(f2) else (x2, f2)
            if fm == 0.0:
                return xm, g.n, True
            tol = 2.0 * 2.22e-16 * abs(xm) + 0.5 * xtol
            span = abs(x2 - x1)
            if span == 0.0:
                return xm, g.n, True
            tl = tol / span
            if tl > 0.5:
                return xm, g.n, True
            inv = 0.0 if math.isinf(m_hat) else 1.0 / m_hat
            G1, G2, G3 = _G(f1, inv), _G(f2, inv), _G(f3, inv)
            t = 0.5
            if x3 != x2 and G3 != G2:
                xi = (x1 - x2) / (x3 - x2)
                ph = (G1 - G2) / (G3 - G2)
                if 0.0 < xi < 1.0 and 1.0 - math.sqrt(1.0 - xi) < ph < math.sqrt(xi):
                    d1, d2, d3 = G1 - G2, G1 - G3, G2 - G3
                    if d1 != 0.0 and d2 != 0.0 and d3 != 0.0:
                        t = (G1 / d1) * (G3 / (-d3)) + \
                            ((x3 - x1) / (x2 - x1)) * (G1 / d2) * (G2 / d3)
            if not math.isfinite(t):
                t = 0.5
            t = min(max(t, tl), 1.0 - tl)
    except _Budget:
        return (x1 if abs(f1) < abs(f2) else x2), g.n, False