#!/usr/bin/env python3
"""
Interference Mapper (R5) -- the FCE wave-engine layer for BSD

Implements the ++/--/+-/-+ interference algebra (FCE reference section 4)
on the AFE phasor walk. Each Euler wave W_p contributes terms
(a_n/n) e^{-2 pi n / sqrt(N)}; the sign of each term against the running
partial sum classifies it:

    ++  term > 0, partial > 0   constructive
    --  term < 0, partial < 0   constructive (troughs reinforce)
    +-  term > 0, partial < 0   destructive
    -+  term < 0, partial > 0   destructive
    0   a_n = 0                 silent prime/index

Rank as null depth -- where the destructive interference REALLY lives:
the completed L-function satisfies Lambda(1) = (1 + eps) * S where
S = sum a_n e^{-x_n}/x_n is the raw phasor walk and eps*S is its FRICKE
BOUNDARY REFLECTION (class-3 event, s <-> 2-s). For eps = -1 the wave
and its reflection cancel EXACTLY: L(E,1) = 0 as a theorem. The raw walk
S itself does NOT vanish (this is precisely why the naive AFE "always
gives nonzero even for rank 1" -- the old pipeline's gotcha). The mapper
therefore detects the null from the certified reflection phase eps, then
switches generator to G_r; the raw-walk null_depth is a descriptive
cancellation metric of the direct wave only.

Null depth := |S_X| / sum |terms| of the raw walk -- near 1 means fully
constructive; markedly lower means heavy internal +-/-+ cancellation
(typical of eps = -1 curves, but the certified null decision is eps).
"""

from dataclasses import dataclass, field
from typing import Dict, List, Optional, Tuple

from flint import arb

from . import certified_arithmetic as ca
from .certified_l_function import CertifiedLFunction, RankCertificate
from .tate_algorithm import analyze_reduction


SIGN_CLASSES = ("++", "--", "+-", "-+", "0")


@dataclass
class PhasorStep:
    n: int
    a_n: int
    term: float                 # (a_n/n) e^{-cn} (midpoint, for profiling)
    partial: float              # running sum midpoint after this term
    sign_class: str             # ++ / -- / +- / -+ / 0


@dataclass
class InterferenceProfile:
    a: int
    b: int
    conductor: int
    epsilon: int                        # certified boundary-reflection phase
    generator_used: int                 # r in G_r after null-switching
    rank_detected: int
    rank_certified: bool
    final_ball: arb                     # certified L^(r)(E,1)/r!
    null_certified: bool                # L(1) ball straddles 0 (or eps=-1 theorem)
    null_depth: float                   # |S_X| / sum|terms| for the r=0 walk
    counts: Dict[str, int]              # sign-class histogram
    constructive_mass: float            # sum |term| over ++/--
    destructive_mass: float             # sum |term| over +-/-+
    steps: List[PhasorStep] = field(default_factory=list, repr=False)

    def fringe(self) -> List[Tuple[int, float]]:
        """(n, partial-sum) sequence -- the fringe profile of the walk."""
        return [(s.n, s.partial) for s in self.steps]


class InterferenceMapper:
    """One FCE instance per Euler wave; composition = the phasor walk."""

    def __init__(self, dps: int = 30):
        self.dps = dps

    def map_curve(self, a: int, b: int,
                  keep_steps: bool = True) -> InterferenceProfile:
        g = analyze_reduction(a, b)
        lf = CertifiedLFunction(a, b, dps=self.dps, reduction=g)

        # --- the r=0 phasor walk (always computed: it IS the fringe) ---
        import math
        N = g.conductor
        X = lf._cutoff()
        an = lf.a_n_list(X)
        c = 2 * math.pi / math.sqrt(N)

        steps: List[PhasorStep] = []
        counts = {k: 0 for k in SIGN_CLASSES}
        partial = 0.0
        cons_mass = 0.0
        dest_mass = 0.0
        for n in range(1, X + 1):
            if an[n] == 0:
                counts["0"] += 1
                continue
            term = 2.0 * an[n] / n * math.exp(-c * n)
            if term > 0:
                cls = "++" if partial >= 0 else "+-"
            else:
                cls = "--" if partial < 0 else "-+"
            partial += term
            counts[cls] += 1
            if cls in ("++", "--"):
                cons_mass += abs(term)
            else:
                dest_mass += abs(term)
            if keep_steps:
                steps.append(PhasorStep(n, an[n], term, partial, cls))

        total_mass = cons_mass + dest_mass
        null_depth = abs(partial) / total_mass if total_mass > 0 else 0.0

        # --- certified null detection + generator switching ---
        cert: RankCertificate = lf.analytic_rank()
        null = cert.rank > 0

        return InterferenceProfile(
            a=a, b=b, conductor=N,
            epsilon=cert.epsilon,
            generator_used=cert.rank,
            rank_detected=cert.rank,
            rank_certified=cert.certified,
            final_ball=cert.leading_coefficient,
            null_certified=null,
            null_depth=null_depth,
            counts=counts,
            constructive_mass=cons_mass,
            destructive_mass=dest_mass,
            steps=steps,
        )

    @staticmethod
    def describe(p: InterferenceProfile) -> str:
        lines = [
            f"Interference map for y^2 = x^3 + {p.a}x + {p.b}  (N={p.conductor})",
            f"  boundary reflection (root number): eps = {p.epsilon:+d}",
            f"  sign-class histogram: " + "  ".join(
                f"{k}:{v}" for k, v in p.counts.items()),
            f"  constructive mass = {p.constructive_mass:.6f}   "
            f"destructive mass = {p.destructive_mass:.6f}",
            f"  null depth |S_X|/sum|terms| = {p.null_depth:.3e}",
        ]
        if p.null_certified:
            lines.append(
                "  DESTRUCTIVE NULL at s=1 (certified): the wave cancels "
                "against its Fricke boundary reflection (eps = -1) -- "
                f"switched generator to G_{p.generator_used}")
        lines.append(
            f"  rank = {p.rank_detected} "
            f"({'certified' if p.rank_certified else 'conditional'}), "
            f"L^({p.rank_detected})(1)/{p.rank_detected}! = {p.final_ball}")
        return "\n".join(lines)