07_phase_transitions.py

"""07 - Phase Transitions: TNFR Bifurcation Dynamics

Exploration of phase transitions and bifurcation behavior in TNFR systems.

PHYSICS: Demonstrates ∂²EPI/∂t² > τ threshold dynamics and controlled bifurcations.
LEARNING: Understanding critical points, hysteresis, and grammar U4 requirements.
"""

import networkx as nx
import numpy as np


def compute_coherence(G):
    """Network phase synchronization: the canonical Kuramoto order
    parameter R = |<e^{iθ}>|.

    R = 1 when phases are fully aligned, R -> 0 when desynchronized
    (random or antiphase). AGENTS.md frames TNFR phase coupling as
    Kuramoto synchronization, so this is the canonical phase-synchrony
    measure. The distinct total coherence
    C(t) = 1/(1 + mean|ΔNFR| + mean|dEPI|) lives in
    tnfr.metrics.coherence and requires the dynamics pipeline.
    """
    thetas = np.array(
        [G.nodes[n].get("theta", G.nodes[n].get("phase", 0.0)) for n in G.nodes()],
        dtype=float,
    )
    if thetas.size == 0:
        return 1.0
    return float(abs(np.mean(np.exp(1j * thetas))))


def compute_delta_nfr(G, node):
    """Compute ΔNFR (structural pressure) for a node."""
    if node not in G.nodes():
        return 0.0

    node_phase = G.nodes[node].get("phase", 0)
    neighbors = list(G.neighbors(node))

    if not neighbors:
        return 0.0

    neighbor_phases = [G.nodes[n].get("phase", 0) for n in neighbors]
    mean_neighbor_phase = np.mean(neighbor_phases)

    phase_diff = abs(node_phase - mean_neighbor_phase)
    return min(phase_diff, 2 * np.pi - phase_diff) / np.pi


def compute_acceleration(G, node, history, dt=0.1):
    """Compute ∂²EPI/∂t² (second derivative) for bifurcation detection."""
    if len(history) < 3:
        return 0.0

    # Use finite difference approximation
    # ∂²φ/∂t² ≈ (φ(t+dt) - 2φ(t) + φ(t-dt)) / dt²

    recent_phases = history[-3:]
    if len(recent_phases) < 3:
        return 0.0

    phi_minus = recent_phases[0]
    phi_current = recent_phases[1]
    phi_plus = recent_phases[2]

    # Handle wraparound for angles
    def angle_diff(a, b):
        diff = a - b
        return np.arctan2(np.sin(diff), np.cos(diff))

    d1 = angle_diff(phi_current, phi_minus) / dt
    d2 = angle_diff(phi_plus, phi_current) / dt

    acceleration = (d2 - d1) / dt
    return abs(acceleration)


def apply_dissonance(G, intensity=0.5):
    """Apply dissonance operator - increases ΔNFR."""
    for node in G.nodes():
        current_phase = G.nodes[node].get("phase", 0)
        # Add random perturbation
        perturbation = np.random.uniform(-intensity, intensity) * np.pi
        G.nodes[node]["phase"] = (current_phase + perturbation) % (2 * np.pi)


def apply_coherence(G, strength=0.3):
    """Apply coherence operator - reduces ΔNFR via negative feedback."""
    # Calculate global phase center
    phases = [G.nodes[n].get("phase", 0) for n in G.nodes()]
    if not phases:
        return

    # Use circular mean for phases
    x = np.mean([np.cos(p) for p in phases])
    y = np.mean([np.sin(p) for p in phases])
    global_phase = np.arctan2(y, x) % (2 * np.pi)

    # Pull all nodes toward global phase
    for node in G.nodes():
        current_phase = G.nodes[node].get("phase", 0)

        # Calculate shortest path to global phase
        diff = global_phase - current_phase
        if diff > np.pi:
            diff -= 2 * np.pi
        elif diff < -np.pi:
            diff += 2 * np.pi

        # Apply coherence correction
        correction = strength * diff
        G.nodes[node]["phase"] = (current_phase + correction) % (2 * np.pi)


def evolve_network_step(G, dt=0.1):
    """Single evolution step with nodal equation."""
    new_phases = {}

    for node in G.nodes():
        current_phase = G.nodes[node].get("phase", 0)
        vf = G.nodes[node].get("vf", 1.0)

        delta_nfr = compute_delta_nfr(G, node)

        neighbors = list(G.neighbors(node))
        if neighbors:
            neighbor_phases = [G.nodes[n].get("phase", 0) for n in neighbors]
            target_phase = np.mean(neighbor_phases)

            direction = target_phase - current_phase
            if direction > np.pi:
                direction -= 2 * np.pi
            elif direction < -np.pi:
                direction += 2 * np.pi

            # Apply nodal equation: ∂EPI/∂t = νf · ΔNFR
            phase_change = vf * delta_nfr * dt * np.sign(direction)
            new_phases[node] = (current_phase + phase_change) % (2 * np.pi)
        else:
            new_phases[node] = current_phase

    for node, phase in new_phases.items():
        G.nodes[node]["phase"] = phase


def phase_transition_experiment():
    """Demonstrate controlled phase transitions and bifurcations."""

    print("🔄 PHASE TRANSITION EXPERIMENT")
    print("━" * 50)
    print("Applying increasing dissonance until bifurcation threshold")

    # Create test network
    G = nx.watts_strogatz_graph(10, 4, 0.1)

    # Initialize with random phases but small νf
    np.random.seed(42)
    for node in G.nodes():
        G.nodes[node]["phase"] = np.random.uniform(0, 2 * np.pi)
        G.nodes[node]["nu_f"] = 0.8
        G.nodes[node]["history"] = []

    initial_coherence = compute_coherence(G)
    print(f"Initial coherence: {initial_coherence:.3f}")

    bifurcation_threshold = 1.5  # τ threshold for ∂²EPI/∂t²
    dissonance_levels = np.linspace(0.1, 1.0, 10)

    results = []

    for level in dissonance_levels:
        print(f"\n🌪️  Dissonance Level: {level:.1f}")

        # Apply dissonance
        apply_dissonance(G, intensity=level)

        # Short evolution to see effect
        phase_history = {node: [] for node in G.nodes()}

        for step in range(15):
            # Record phases before step
            for node in G.nodes():
                phase_history[node].append(G.nodes[node]["phase"])

            evolve_network_step(G)

        # Check for bifurcations
        bifurcations_detected = 0
        max_acceleration = 0

        for node in G.nodes():
            history = phase_history[node]
            acceleration = compute_acceleration(G, node, history)
            max_acceleration = max(max_acceleration, acceleration)

            if acceleration > bifurcation_threshold:
                bifurcations_detected += 1

        final_coherence = compute_coherence(G)
        avg_delta_nfr = np.mean([compute_delta_nfr(G, n) for n in G.nodes()])

        print(f"   Final coherence: {final_coherence:.3f}")
        print(f"   Max acceleration: {max_acceleration:.3f}")
        print(f"   Bifurcations detected: {bifurcations_detected}")
        print(f"   Average ΔNFR: {avg_delta_nfr:.3f}")

        results.append(
            {
                "level": level,
                "coherence": final_coherence,
                "acceleration": max_acceleration,
                "bifurcations": bifurcations_detected,
                "delta_nfr": avg_delta_nfr,
            }
        )

        # Apply coherence to stabilize (Grammar U4a requirement)
        if bifurcations_detected > 0:
            print("   🛡️  Applying coherence (U4a stabilizer)")
            apply_coherence(G, strength=0.5)
            stabilized_coherence = compute_coherence(G)
            print(f"   Stabilized coherence: {stabilized_coherence:.3f}")

    return results


def hysteresis_experiment():
    """Demonstrate hysteresis in phase transitions."""

    print("\n🔄 HYSTERESIS EXPERIMENT")
    print("━" * 50)
    print("Forward: Increasing dissonance")
    print("Backward: Decreasing dissonance")

    G = nx.cycle_graph(8)

    # Initialize synchronized state
    for node in G.nodes():
        G.nodes[node]["phase"] = 0.0  # All synchronized
        G.nodes[node]["nu_f"] = 1.0

    # Forward path: increasing dissonance
    dissonance_levels = np.linspace(0.0, 0.8, 16)
    forward_coherence = []

    for level in dissonance_levels:
        apply_dissonance(G, intensity=level)

        # Brief evolution
        for _ in range(10):
            evolve_network_step(G)

        coherence = compute_coherence(G)
        forward_coherence.append(coherence)

    # Backward path: decreasing dissonance
    backward_coherence = []

    for level in reversed(dissonance_levels):
        # Apply inverse of dissonance (coherence)
        apply_coherence(G, strength=0.8 - level)

        # Brief evolution
        for _ in range(10):
            evolve_network_step(G)

        coherence = compute_coherence(G)
        backward_coherence.append(coherence)

    # Analyze hysteresis
    print("Hysteresis Analysis:")
    for i in range(len(dissonance_levels)):
        level = dissonance_levels[i]
        forward_c = forward_coherence[i]
        backward_c = backward_coherence[-(i + 1)]  # Reverse order
        hysteresis = abs(forward_c - backward_c)

        if i % 4 == 0:  # Print every 4th point
            print(
                f"  Level {level:.2f}: Forward {forward_c:.3f}, "
                f"Backward {backward_c:.3f}, Δ = {hysteresis:.3f}"
            )

    avg_hysteresis = np.mean(
        [
            abs(forward_coherence[i] - backward_coherence[-(i + 1)])
            for i in range(len(dissonance_levels))
        ]
    )

    print(f"Average hysteresis: {avg_hysteresis:.3f}")

    if avg_hysteresis > 0.05:
        print("🔄 Strong hysteresis detected - system exhibits memory")
    elif avg_hysteresis > 0.02:
        print("🔄 Moderate hysteresis - weak memory effects")
    else:
        print("🔄 Minimal hysteresis - system is reversible")


def critical_point_analysis():
    """Find critical points in parameter space."""

    print("\n🎯 CRITICAL POINT ANALYSIS")
    print("━" * 50)

    # Test different network sizes
    sizes = [6, 8, 10, 12, 16]
    critical_points = []

    for size in sizes:
        print(f"\nNetwork size: {size}")

        G = nx.cycle_graph(size)

        # Initialize random state
        np.random.seed(42 + size)
        for node in G.nodes():
            G.nodes[node]["phase"] = np.random.uniform(0, 2 * np.pi)
            G.nodes[node]["nu_f"] = 1.0

        # Find critical dissonance level
        critical_level = None

        for level in np.linspace(0.1, 1.0, 20):
            # Reset network
            for node in G.nodes():
                G.nodes[node]["phase"] = np.random.uniform(0, 2 * np.pi)

            apply_dissonance(G, intensity=level)

            # Measure response
            initial_coherence = compute_coherence(G)

            # Short evolution
            for _ in range(20):
                evolve_network_step(G)

            final_coherence = compute_coherence(G)
            response = abs(final_coherence - initial_coherence)

            # Look for sharp response increase (critical point)
            if response > 0.3:  # Threshold for significant response
                critical_level = level
                break

        if critical_level:
            print(f"  Critical level: ~{critical_level:.2f}")
            critical_points.append(critical_level)
        else:
            print("  No clear critical point found")

    if critical_points:
        avg_critical = np.mean(critical_points)
        std_critical = np.std(critical_points)
        print(f"\nCritical point statistics:")
        print(f"  Average: {avg_critical:.3f} ± {std_critical:.3f}")
        print(f"  Range: [{min(critical_points):.3f}, {max(critical_points):.3f}]")


def grammar_u4_demonstration():
    """Demonstrate Grammar U4 requirements for bifurcation handling."""

    print("\n📋 GRAMMAR U4 DEMONSTRATION")
    print("━" * 50)
    print("U4a: Bifurcation triggers need handlers")
    print("U4b: Transformers need context")

    G = nx.complete_graph(6)

    # Initialize synchronized
    for node in G.nodes():
        G.nodes[node]["phase"] = 0.0
        G.nodes[node]["nu_f"] = 1.0

    print("\n1️⃣ INCORRECT: Dissonance without handler")

    initial_coherence = compute_coherence(G)
    print(f"   Initial coherence: {initial_coherence:.3f}")

    # Apply strong dissonance (violates U4a)
    apply_dissonance(G, intensity=0.8)

    # Let it evolve without stabilizer
    for _ in range(20):
        evolve_network_step(G)

    uncontrolled_coherence = compute_coherence(G)
    print(f"   Uncontrolled result: {uncontrolled_coherence:.3f}")
    print("   ❌ Violates U4a - no handler for bifurcation trigger")

    print("\n2️⃣ CORRECT: Dissonance with coherence handler")

    # Reset
    for node in G.nodes():
        G.nodes[node]["phase"] = 0.0

    initial_coherence = compute_coherence(G)
    print(f"   Initial coherence: {initial_coherence:.3f}")

    # Apply dissonance (trigger)
    apply_dissonance(G, intensity=0.8)

    # Brief evolution
    for _ in range(10):
        evolve_network_step(G)

    mid_coherence = compute_coherence(G)
    print(f"   After dissonance: {mid_coherence:.3f}")

    # Apply coherence (handler - satisfies U4a)
    apply_coherence(G, strength=0.6)

    # Continue evolution
    for _ in range(10):
        evolve_network_step(G)

    controlled_coherence = compute_coherence(G)
    print(f"   After coherence: {controlled_coherence:.3f}")
    print("   ✅ Satisfies U4a - handler controls bifurcation")

    improvement = controlled_coherence - uncontrolled_coherence
    print(f"   Improvement: {improvement:+.3f}")

    print("\n🔬 U4 Physics Insight:")
    print("   • Dissonance creates ∂²EPI/∂t² > τ")
    print("   • Without handlers → chaos/fragmentation")
    print("   • With handlers → controlled reorganization")
    print("   • Grammar U4 = physics constraint, not arbitrary rule")


def phase_transitions_demo():
    """Comprehensive demonstration of TNFR phase transitions."""

    print("=" * 80)
    print("                🔄 PHASE TRANSITION DYNAMICS 🔄")
    print("=" * 80)
    print()
    print("Exploring bifurcations, critical points, and Grammar U4 physics...")
    print("PHYSICS: ∂²EPI/∂t² > τ triggers require stabilizers (U4a)")
    print("INSIGHT: Phase transitions reveal deep structure-dynamics coupling")
    print()

    # Run experiments
    transition_results = phase_transition_experiment()

    hysteresis_experiment()

    critical_point_analysis()

    grammar_u4_demonstration()

    # Summary analysis
    print("\n" + "=" * 80)
    print("🧮 PHASE TRANSITION INSIGHTS")
    print("=" * 80)

    # Find transition point
    coherence_drop_threshold = 0.2
    transition_point = None

    for result in transition_results:
        if result["coherence"] < (
            transition_results[0]["coherence"] - coherence_drop_threshold
        ):
            transition_point = result["level"]
            break

    if transition_point:
        print(f"📊 Coherence transition at dissonance level: {transition_point:.1f}")
    else:
        print("📊 No sharp coherence transition detected")

    # Bifurcation statistics
    total_bifurcations = sum(r["bifurcations"] for r in transition_results)
    max_bifurcations = max(r["bifurcations"] for r in transition_results)

    print(f"⚡ Total bifurcations detected: {total_bifurcations}")
    print(f"⚡ Maximum bifurcations in single experiment: {max_bifurcations}")

    # Universal principles
    print("\n🔬 UNIVERSAL PHASE TRANSITION PRINCIPLES:")
    print("━" * 60)
    print("• Dissonance increases ∂²EPI/∂t² (acceleration)")
    print("• Threshold τ triggers bifurcation dynamics")
    print("• Uncontrolled bifurcations → fragmentation")
    print("• Stabilizers (coherence) → controlled reorganization")
    print("• Hysteresis → system memory and path dependence")
    print("• Critical points → sharp response transitions")
    print("• Network size affects critical thresholds")

    print("\n🛡️ GRAMMAR U4 VALIDATION:")
    print("━" * 60)
    print("✅ U4a enforced: Bifurcation triggers paired with handlers")
    print("✅ U4b respected: Transformers require destabilizer context")
    print("✅ Physics basis: ∫νf·ΔNFR dt convergence requires stabilizers")
    print("✅ Practical result: Controlled vs chaotic reorganization")

    print("\n🚀 ADVANCED PHASE TRANSITION RESEARCH:")
    print("━" * 60)
    print("• Multi-parameter phase diagrams (νf, coupling strength)")
    print("• Avalanche dynamics and self-organized criticality")
    print("• Phase transitions in hierarchical networks")
    print("• Temperature-like parameters for thermal transitions")
    print("• Quantum-inspired coherent/decoherent phases")


if __name__ == "__main__":
    phase_transitions_demo()