test_variational.py

"""Tests for TNFR Variational Principle — Lagrangian Action Formulation.

Validates that the nodal equation ∂EPI/∂t = νf · ΔNFR(t) arises from
the Euler-Lagrange equations of the TNFR action functional:

    S_TNFR = ∫ dt Σ_i ℒ_TNFR(i)

where ℒ = T − V with T = ½(J_φ² + J_ΔNFR²) and V = ½(Φ_s² + |∇φ|² + K_φ²).

Tests verify:
1.  Lagrangian density: ℒ = T − V (sign and magnitude)
2.  Hamiltonian density: H = T + V = ½ · energy_density (consistency)
3.  Conjugate pairs: (K_φ, J_φ) and (Φ_s, J_ΔNFR)
4.  Euler-Lagrange residual: small for grammar-compliant evolution
5.  Action functional: finite for U2-compliant sequences
6.  Symplectic preservation: canonical operators preserve ω
7.  Grammar as stationarity: U1-U6 mapped to variational conditions
8.  Potential critical points: thresholds at φ, γ/π, 0.9π
9.  VariationalTracker: time-series accumulation
10. Operator classification: generating/dissipative/canonical
11. Cross-topology validation: WS, BA, Grid
12. Consistency with conservation.py energy functional
13. Virial ratio diagnostics
14. Reproducibility under deterministic seeds

TIER: CORE PHYSICS — variational formulation axiomatises the nodal equation.
"""

from __future__ import annotations

import copy
import math
import os
import sys

import networkx as nx
import numpy as np
import pytest

sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "..", "src"))

from tnfr.constants import inject_defaults
from tnfr.constants.canonical import PI, U6_STRUCTURAL_POTENTIAL_LIMIT
from tnfr.physics.conservation import compute_energy_functional
from tnfr.physics.unified import compute_energy_density
from tnfr.physics.variational import (
    ConjugatePair,
    CriticalPointAnalysis,
    EulerLagrangeResidual,
    GrammarStationarityAnalysis,
    LagrangianSnapshot,
    SymplecticCheck,
    VariationalTimeSeries,
    VariationalTracker,
    analyze_grammar_stationarity,
    analyze_potential_critical_points,
    capture_lagrangian_snapshot,
    check_symplectic_preservation,
    classify_operator_canonical,
    compute_action_functional,
    compute_euler_lagrange_residual,
    compute_hamiltonian_density,
    compute_interaction_density,
    compute_kinetic_density,
    compute_lagrangian_density,
    compute_phase_space_volume,
    compute_poisson_bracket_estimate,
    compute_potential_density,
    compute_variational_suite,
    identify_conjugate_pairs,
)

# ---------------------------------------------------------------------------
# Fixtures
# ---------------------------------------------------------------------------


def _make_tnfr_graph(
    n: int = 30,
    topology: str = "watts_strogatz",
    seed: int = 42,
) -> nx.Graph:
    """Build a TNFR-ready graph with canonical attributes."""
    rng = np.random.default_rng(seed)

    if topology == "watts_strogatz":
        G = nx.watts_strogatz_graph(n, 4, 0.3, seed=seed)
    elif topology == "barabasi_albert":
        G = nx.barabasi_albert_graph(n, 3, seed=seed)
    elif topology == "grid":
        side = int(math.sqrt(n))
        G = nx.grid_2d_graph(side, side)
    else:
        G = nx.watts_strogatz_graph(n, 4, 0.3, seed=seed)

    inject_defaults(G)

    for node in G.nodes():
        G.nodes[node]["phase"] = rng.uniform(0, 2 * math.pi)
        G.nodes[node]["frequency"] = rng.uniform(0.1, 1.0)
        G.nodes[node]["delta_nfr"] = rng.uniform(-0.5, 0.5)
        G.nodes[node]["EPI"] = f"epi_{node}"

    return G


def _perturb_graph(G: nx.Graph, seed: int = 99) -> nx.Graph:
    """Create a slightly perturbed copy for two-snapshot tests."""
    G2 = copy.deepcopy(G)
    rng = np.random.default_rng(seed)
    for node in G2.nodes():
        G2.nodes[node]["phase"] += rng.uniform(-0.1, 0.1)
        G2.nodes[node]["delta_nfr"] += rng.uniform(-0.05, 0.05)
    return G2


@pytest.fixture
def ws_graph():
    return _make_tnfr_graph(30, "watts_strogatz", seed=42)


@pytest.fixture
def ba_graph():
    return _make_tnfr_graph(30, "barabasi_albert", seed=42)


@pytest.fixture
def grid_graph():
    return _make_tnfr_graph(25, "grid", seed=42)


# ---------------------------------------------------------------------------
# 1. Lagrangian density: ℒ = T − V
# ---------------------------------------------------------------------------


class TestLagrangianDensity:
    """ℒ(i) = T(i) − V(i) with correct sign and magnitude."""

    def test_lagrangian_equals_T_minus_V(self, ws_graph):
        T = compute_kinetic_density(ws_graph)
        V = compute_potential_density(ws_graph)
        L = compute_lagrangian_density(ws_graph)
        for n in ws_graph.nodes():
            assert abs(L[n] - (T[n] - V[n])) < 1e-12

    def test_T_non_negative(self, ws_graph):
        T = compute_kinetic_density(ws_graph)
        for v in T.values():
            assert v >= 0.0

    def test_V_non_negative(self, ws_graph):
        V = compute_potential_density(ws_graph)
        for v in V.values():
            assert v >= 0.0

    def test_lagrangian_can_be_negative(self, ws_graph):
        """Potential-dominated states have ℒ < 0 (attractor basins)."""
        L = compute_lagrangian_density(ws_graph)
        # At least some nodes should have negative Lagrangian
        vals = list(L.values())
        assert any(v < 0 for v in vals) or any(v > 0 for v in vals)


# ---------------------------------------------------------------------------
# 2. Hamiltonian density: H = T + V = ½ · energy_density
# ---------------------------------------------------------------------------


class TestHamiltonianConsistency:
    """H(i) = T(i) + V(i) and H = ½ℰ from unified.py."""

    def test_hamiltonian_equals_T_plus_V(self, ws_graph):
        T = compute_kinetic_density(ws_graph)
        V = compute_potential_density(ws_graph)
        H = compute_hamiltonian_density(ws_graph)
        for n in ws_graph.nodes():
            assert abs(H[n] - (T[n] + V[n])) < 1e-12

    def test_hamiltonian_half_energy_density(self, ws_graph):
        """H(i) = ½ · ℰ(i) where ℰ is from unified.compute_energy_density."""
        H = compute_hamiltonian_density(ws_graph)
        E = compute_energy_density(ws_graph)
        for n in ws_graph.nodes():
            assert (
                abs(H[n] - 0.5 * E[n]) < 1e-12
            ), f"Node {n}: H={H[n]:.6f}, ½ℰ={0.5*E[n]:.6f}"

    def test_total_hamiltonian_equals_energy_functional(self, ws_graph):
        """Σ H(i) = compute_energy_functional(G)."""
        H = compute_hamiltonian_density(ws_graph)
        total_H = sum(H.values())
        E_func = compute_energy_functional(ws_graph)
        assert (
            abs(total_H - E_func) < 1e-10
        ), f"Total H={total_H:.6f}, E_func={E_func:.6f}"


# ---------------------------------------------------------------------------
# 3. Conjugate pairs
# ---------------------------------------------------------------------------


class TestConjugatePairs:
    """Canonical conjugate pairs: (K_φ, J_φ) and (Φ_s, J_ΔNFR)."""

    def test_two_sectors_identified(self, ws_graph):
        geo, pot = identify_conjugate_pairs(ws_graph)
        assert geo.sector == "geometric"
        assert pot.sector == "potential"

    def test_pairs_have_matching_nodes(self, ws_graph):
        geo, pot = identify_conjugate_pairs(ws_graph)
        nodes = set(ws_graph.nodes())
        assert set(geo.q.keys()) == nodes
        assert set(geo.p.keys()) == nodes
        assert set(pot.q.keys()) == nodes
        assert set(pot.p.keys()) == nodes

    def test_poisson_bracket_non_degenerate(self, ws_graph):
        """Non-zero Poisson bracket indicates non-degenerate symplectic structure."""
        geo, pot = identify_conjugate_pairs(ws_graph)
        pb_geo = compute_poisson_bracket_estimate(geo)
        pb_pot = compute_poisson_bracket_estimate(pot)
        # For a random graph with varied fields, brackets should be non-zero
        assert pb_geo > 0.0 or pb_pot > 0.0

    def test_phase_space_volume_positive(self, ws_graph):
        geo, pot = identify_conjugate_pairs(ws_graph)
        vol_geo = compute_phase_space_volume(geo)
        vol_pot = compute_phase_space_volume(pot)
        assert vol_geo >= 0.0
        assert vol_pot >= 0.0


# ---------------------------------------------------------------------------
# 4. Euler-Lagrange residual
# ---------------------------------------------------------------------------


class TestEulerLagrangeResidual:
    """EL residual quantifies departure from stationarity."""

    def test_equilibrium_has_small_residual(self, ws_graph):
        """Identical snapshots → zero residual."""
        snap = capture_lagrangian_snapshot(ws_graph)
        el = compute_euler_lagrange_residual(snap, snap, dt=1.0)
        # Identical snapshots: dp/dt = 0, but avg q remains → residual = |q_avg|
        # Still, it should be finite and well-defined
        assert math.isfinite(el.rms_residual)
        assert el.stationarity_quality > 0

    def test_perturbed_has_larger_residual(self, ws_graph):
        """Perturbation increases EL residual."""
        snap_before = capture_lagrangian_snapshot(ws_graph)
        G2 = _perturb_graph(ws_graph, seed=99)
        snap_after = capture_lagrangian_snapshot(G2)
        el = compute_euler_lagrange_residual(snap_before, snap_after, dt=1.0)
        assert el.rms_residual >= 0.0
        assert 0.0 < el.stationarity_quality <= 1.0

    def test_residual_structure(self, ws_graph):
        snap = capture_lagrangian_snapshot(ws_graph)
        G2 = _perturb_graph(ws_graph)
        snap2 = capture_lagrangian_snapshot(G2)
        el = compute_euler_lagrange_residual(snap, snap2)
        assert isinstance(el, EulerLagrangeResidual)
        assert len(el.residual) == ws_graph.number_of_nodes()
        assert el.max_residual >= el.rms_residual >= el.mean_residual >= 0


# ---------------------------------------------------------------------------
# 5. Action functional
# ---------------------------------------------------------------------------


class TestActionFunctional:
    """S = ∫ dt L is finite for well-behaved sequences."""

    def test_single_snapshot_zero_action(self, ws_graph):
        snap = capture_lagrangian_snapshot(ws_graph)
        S = compute_action_functional([snap], dt=1.0)
        assert math.isfinite(S)

    def test_two_snapshots_finite_action(self, ws_graph):
        snap1 = capture_lagrangian_snapshot(ws_graph)
        G2 = _perturb_graph(ws_graph)
        snap2 = capture_lagrangian_snapshot(G2)
        S = compute_action_functional([snap1, snap2], dt=1.0)
        assert math.isfinite(S)

    def test_action_scales_with_dt(self, ws_graph):
        snap = capture_lagrangian_snapshot(ws_graph)
        S1 = compute_action_functional([snap], dt=1.0)
        S2 = compute_action_functional([snap], dt=2.0)
        assert abs(S2 - 2.0 * S1) < 1e-12


# ---------------------------------------------------------------------------
# 6. Symplectic preservation
# ---------------------------------------------------------------------------


class TestSymplecticPreservation:
    """Canonical operators preserve symplectic structure."""

    def test_identity_is_canonical(self, ws_graph):
        """No change → perfect canonical transformation."""
        snap = capture_lagrangian_snapshot(ws_graph)
        sc = check_symplectic_preservation(snap, snap, "identity")
        assert sc.is_canonical
        assert sc.classification == "canonical"
        assert abs(sc.volume_ratio - 1.0) < 1e-12

    def test_perturbation_classified(self, ws_graph):
        snap1 = capture_lagrangian_snapshot(ws_graph)
        G2 = _perturb_graph(ws_graph)
        snap2 = capture_lagrangian_snapshot(G2)
        sc = check_symplectic_preservation(snap1, snap2, "perturbation")
        assert isinstance(sc, SymplecticCheck)
        assert sc.classification in ("canonical", "dissipative", "expansive", "mixed")

    def test_symplectic_check_structure(self, ws_graph):
        snap = capture_lagrangian_snapshot(ws_graph)
        sc = check_symplectic_preservation(snap, snap, "test")
        assert sc.operator_name == "test"
        assert math.isfinite(sc.symplectic_ratio_geometric)
        assert math.isfinite(sc.symplectic_ratio_potential)


# ---------------------------------------------------------------------------
# 7. Grammar as stationarity conditions
# ---------------------------------------------------------------------------


class TestGrammarStationarity:
    """Grammar rules U1-U6 mapped to variational conditions."""

    def test_all_six_rules_covered(self, ws_graph):
        results = analyze_grammar_stationarity(ws_graph)
        rules = {r.rule for r in results}
        assert "U1a" in rules
        assert "U1b" in rules
        assert "U2" in rules
        assert "U3" in rules
        assert "U4" in rules
        assert "U5" in rules
        assert "U6" in rules

    def test_each_has_interpretation(self, ws_graph):
        results = analyze_grammar_stationarity(ws_graph)
        for r in results:
            assert len(r.variational_interpretation) > 10
            assert isinstance(r.is_satisfied, bool)
            assert math.isfinite(r.diagnostic_value)

    def test_with_snapshots(self, ws_graph):
        snap1 = capture_lagrangian_snapshot(ws_graph)
        G2 = _perturb_graph(ws_graph)
        snap2 = capture_lagrangian_snapshot(G2)
        results = analyze_grammar_stationarity(G2, snapshots=[snap1, snap2], dt=1.0)
        # U2 should use action-based check when snapshots provided
        u2 = [r for r in results if r.rule == "U2"][0]
        assert math.isfinite(u2.diagnostic_value)


# ---------------------------------------------------------------------------
# 8. Potential critical points (thresholds)
# ---------------------------------------------------------------------------


class TestCriticalPoints:
    """TNFR thresholds correspond to critical points of V."""

    def test_three_fields_analysed(self, ws_graph):
        results = analyze_potential_critical_points(ws_graph)
        names = {r.field_name for r in results}
        assert "Phi_s" in names
        assert "grad_phi" in names
        assert "K_phi" in names

    def test_thresholds_match_theory(self, ws_graph):
        results = analyze_potential_critical_points(ws_graph)
        for r in results:
            if r.field_name == "Phi_s":
                assert abs(r.threshold_value - U6_STRUCTURAL_POTENTIAL_LIMIT) < 1e-10
            elif r.field_name == "grad_phi":
                assert abs(r.threshold_value - 0.9 * PI) < 1e-10
            elif r.field_name == "K_phi":
                assert abs(r.threshold_value - 0.9 * PI) < 1e-10

    def test_critical_type_valid(self, ws_graph):
        results = analyze_potential_critical_points(ws_graph)
        for r in results:
            assert r.critical_type in ("minimum", "maximum", "saddle", "regular")


# ---------------------------------------------------------------------------
# 9. VariationalTracker
# ---------------------------------------------------------------------------


class TestVariationalTracker:
    """Time-series tracker for variational diagnostics."""

    def test_single_record(self, ws_graph):
        tracker = VariationalTracker(ws_graph)
        snap = tracker.record(t=0.0)
        assert isinstance(snap, LagrangianSnapshot)
        report = tracker.report()
        assert len(report.times) == 1
        assert report.total_lagrangian[0] == snap.total_lagrangian

    def test_two_records(self, ws_graph):
        tracker = VariationalTracker(ws_graph)
        tracker.record(t=0.0)
        # Perturb
        for n in ws_graph.nodes():
            ws_graph.nodes[n]["phase"] += 0.01
        tracker.record(t=1.0)
        report = tracker.report()
        assert len(report.times) == 2
        assert len(report.el_rms_residual) == 2
        assert report.el_rms_residual[0] == 0.0  # first step
        assert report.el_rms_residual[1] >= 0.0  # second step

    def test_action_accumulation(self, ws_graph):
        tracker = VariationalTracker(ws_graph)
        tracker.record(t=0.0)
        for n in ws_graph.nodes():
            ws_graph.nodes[n]["delta_nfr"] *= 0.9
        tracker.record(t=1.0)
        assert math.isfinite(tracker.action)
        report = tracker.report()
        assert report.is_action_finite

    def test_latest_snapshot(self, ws_graph):
        tracker = VariationalTracker(ws_graph)
        assert tracker.latest_snapshot is None
        snap = tracker.record(t=0.0)
        assert tracker.latest_snapshot is snap


# ---------------------------------------------------------------------------
# 10. Operator classification
# ---------------------------------------------------------------------------


class TestOperatorClassification:
    """Classify operators as generating/dissipative/canonical."""

    def test_identity_classified_neutral(self, ws_graph):
        snap = capture_lagrangian_snapshot(ws_graph)
        result = classify_operator_canonical(snap, snap, "SHA")
        assert result["energy_classification"] == "neutral"
        assert result["consistent_with_theory"]

    def test_energy_increase_classified_generating(self, ws_graph):
        snap_before = capture_lagrangian_snapshot(ws_graph)
        # Increase energy by amplifying ΔNFR
        G2 = copy.deepcopy(ws_graph)
        for n in G2.nodes():
            G2.nodes[n]["delta_nfr"] *= 3.0
        snap_after = capture_lagrangian_snapshot(G2)
        result = classify_operator_canonical(snap_before, snap_after, "OZ")
        # OZ should increase energy
        assert (
            result["energy_change"] > 0 or result["energy_classification"] == "neutral"
        )

    def test_all_operators_have_expected_type(self):
        """Every canonical operator has a theoretical classification."""
        from tnfr.physics.variational import _OPERATOR_CANONICAL_MAP

        expected_ops = {
            "AL",
            "EN",
            "IL",
            "OZ",
            "UM",
            "RA",
            "SHA",
            "VAL",
            "NUL",
            "THOL",
            "ZHIR",
            "NAV",
            "REMESH",
        }
        assert set(_OPERATOR_CANONICAL_MAP.keys()) == expected_ops


# ---------------------------------------------------------------------------
# 11. Cross-topology validation
# ---------------------------------------------------------------------------


class TestCrossTopology:
    """Variational principle holds across topologies."""

    @pytest.mark.parametrize(
        "topology,n",
        [
            ("watts_strogatz", 30),
            ("barabasi_albert", 30),
            ("grid", 25),
        ],
    )
    def test_lagrangian_defined(self, topology, n):
        G = _make_tnfr_graph(n, topology, seed=42)
        L = compute_lagrangian_density(G)
        assert len(L) == G.number_of_nodes()
        assert all(math.isfinite(v) for v in L.values())

    @pytest.mark.parametrize(
        "topology,n",
        [
            ("watts_strogatz", 30),
            ("barabasi_albert", 30),
            ("grid", 25),
        ],
    )
    def test_hamiltonian_consistent(self, topology, n):
        G = _make_tnfr_graph(n, topology, seed=42)
        H = compute_hamiltonian_density(G)
        E_func = compute_energy_functional(G)
        assert abs(sum(H.values()) - E_func) < 1e-10

    @pytest.mark.parametrize(
        "topology,n",
        [
            ("watts_strogatz", 30),
            ("barabasi_albert", 30),
            ("grid", 25),
        ],
    )
    def test_variational_suite(self, topology, n):
        G = _make_tnfr_graph(n, topology, seed=42)
        suite = compute_variational_suite(G)
        assert "lagrangian_snapshot" in suite
        assert "critical_points" in suite
        assert "grammar_stationarity" in suite
        assert math.isfinite(suite["virial_ratio"])


# ---------------------------------------------------------------------------
# 12. Consistency with conservation.py
# ---------------------------------------------------------------------------


class TestConservationConsistency:
    """Variational module consistent with conservation module."""

    def test_energy_functional_consistent(self, ws_graph):
        """Total Hamiltonian = conservation energy functional."""
        snap = capture_lagrangian_snapshot(ws_graph)
        E_cons = compute_energy_functional(ws_graph)
        assert abs(snap.total_hamiltonian - E_cons) < 1e-10

    def test_conjugate_pairs_match_conservation_sectors(self, ws_graph):
        """Geometric/potential sectors match conservation decomposition."""
        snap = capture_lagrangian_snapshot(ws_graph)
        # Geometric: q = K_φ, p = J_φ → conservation geometric sector
        # Potential: q = Φ_s, p = J_ΔNFR → conservation potential sector
        from tnfr.physics.conservation import compute_charge_density

        rho = compute_charge_density(ws_graph)
        # ρ = Φ_s + K_φ
        for n in ws_graph.nodes():
            rho_from_pairs = (
                snap.conjugate_potential.q[n] + snap.conjugate_geometric.q[n]
            )
            assert abs(rho[n] - rho_from_pairs) < 1e-12


# ---------------------------------------------------------------------------
# 13. Virial ratio and energy partition
# ---------------------------------------------------------------------------


class TestVirialRatio:
    """Virial ratio T/V diagnostics."""

    def test_virial_computable(self, ws_graph):
        suite = compute_variational_suite(ws_graph)
        assert math.isfinite(suite["virial_ratio"])
        assert suite["virial_ratio"] >= 0


# ---------------------------------------------------------------------------
# 14. Reproducibility
# ---------------------------------------------------------------------------


class TestReproducibility:
    """Deterministic seeds → identical results."""

    def test_same_seed_same_lagrangian(self):
        G1 = _make_tnfr_graph(30, "watts_strogatz", seed=42)
        G2 = _make_tnfr_graph(30, "watts_strogatz", seed=42)
        L1 = compute_lagrangian_density(G1)
        L2 = compute_lagrangian_density(G2)
        for n in G1.nodes():
            assert abs(L1[n] - L2[n]) < 1e-14

    def test_same_seed_same_action(self):
        G1 = _make_tnfr_graph(30, "watts_strogatz", seed=42)
        G2 = _make_tnfr_graph(30, "watts_strogatz", seed=42)
        snap1 = capture_lagrangian_snapshot(G1)
        snap2 = capture_lagrangian_snapshot(G2)
        S1 = compute_action_functional([snap1], dt=1.0)
        S2 = compute_action_functional([snap2], dt=1.0)
        assert abs(S1 - S2) < 1e-14


# ---------------------------------------------------------------------------
# 15. Interaction density (bilinear coupling)
# ---------------------------------------------------------------------------


class TestInteractionDensity:
    """Interaction ℒ_int = existing action_density."""

    def test_interaction_matches_action_density(self, ws_graph):
        from tnfr.physics.unified import compute_action_density

        interaction = compute_interaction_density(ws_graph)
        action_d = compute_action_density(ws_graph)
        for n in ws_graph.nodes():
            assert abs(interaction[n] - action_d[n]) < 1e-12


# ---------------------------------------------------------------------------
# 16. Snapshot completeness
# ---------------------------------------------------------------------------


class TestSnapshotCompleteness:
    """LagrangianSnapshot contains all required information."""

    def test_snapshot_fields(self, ws_graph):
        snap = capture_lagrangian_snapshot(ws_graph)
        assert isinstance(snap, LagrangianSnapshot)
        N = ws_graph.number_of_nodes()
        assert len(snap.kinetic) == N
        assert len(snap.potential) == N
        assert len(snap.lagrangian) == N
        assert len(snap.hamiltonian) == N
        assert len(snap.interaction) == N
        assert math.isfinite(snap.total_lagrangian)
        assert math.isfinite(snap.total_hamiltonian)
        assert snap.total_hamiltonian >= 0  # energy ≥ 0

    def test_snapshot_totals_consistent(self, ws_graph):
        snap = capture_lagrangian_snapshot(ws_graph)
        assert (
            abs(snap.total_lagrangian - (snap.total_kinetic - snap.total_potential))
            < 1e-12
        )
        assert (
            abs(snap.total_hamiltonian - (snap.total_kinetic + snap.total_potential))
            < 1e-12
        )


# ---------------------------------------------------------------------------
# 17. Sector translation: variational ↔ conservation ↔ unified
# ---------------------------------------------------------------------------


class TestSectorTranslation:
    """translate_sectors() bridges the three decompositions of the 6D field."""

    def test_keys_present(self, ws_graph):
        from tnfr.physics.variational import translate_sectors

        result = translate_sectors(ws_graph)
        assert "variational" in result
        assert "conservation" in result
        assert "unified_psi" in result
        assert "energy_density" in result
        assert "consistency_check" in result

    def test_consistency_check_near_zero(self, ws_graph):
        """T(i) + V(i) == ½·ℰ(i) for every node."""
        from tnfr.physics.variational import translate_sectors

        result = translate_sectors(ws_graph)
        assert result["consistency_check"] < 1e-12

    def test_variational_sector_node_coverage(self, ws_graph):
        from tnfr.physics.variational import translate_sectors

        result = translate_sectors(ws_graph)
        nodes = set(ws_graph.nodes())
        assert set(result["variational"]["T"].keys()) == nodes
        assert set(result["variational"]["V"].keys()) == nodes

    def test_conservation_sector_node_coverage(self, ws_graph):
        from tnfr.physics.variational import translate_sectors

        result = translate_sectors(ws_graph)
        nodes = set(ws_graph.nodes())
        assert set(result["conservation"]["rho"].keys()) == nodes
        assert set(result["conservation"]["J_phi"].keys()) == nodes
        assert set(result["conservation"]["J_dnfr"].keys()) == nodes

    def test_unified_psi_is_complex(self, ws_graph):
        from tnfr.physics.variational import translate_sectors

        result = translate_sectors(ws_graph)
        for psi in result["unified_psi"].values():
            assert isinstance(psi, complex)

    def test_energy_density_matches_unified(self, ws_graph):
        """Raw ℰ returned by translate_sectors matches unified.py directly."""
        from tnfr.physics.variational import translate_sectors

        result = translate_sectors(ws_graph)
        raw_direct = compute_energy_density(ws_graph)
        for n in ws_graph.nodes():
            assert abs(result["energy_density"][n] - raw_direct[n]) < 1e-12

    def test_cross_topology_consistency(self, ba_graph, grid_graph):
        """Sector translation holds across BA and grid topologies."""
        from tnfr.physics.variational import translate_sectors

        for G in (ba_graph, grid_graph):
            result = translate_sectors(G)
            assert result["consistency_check"] < 1e-12