#!/usr/bin/env python3 """ B.6.7 SELF-CONSISTENT f+A|t UNFOLDING TEST ============================================= Date: 2026-03-19 Purpose: Test whether a pulsed N=3 FDTD simulation with SELF-CONSISTENT feedback (accumulated intensity → potential → decoherence ratio) produces emergent curvature, pentagonal defects, and chirality. THEORETICAL BASIS: 1. f|t produces hexagonal interference (N=3, proven A.1) 2. C_potential breaks symmetry via decoherence (proven, AUDITED) 3. Energy coalescence deepens the potential (theory.txt lines 43-46) 4. The potential is SCALAR — it grows with accumulated energy 5. The geometry IS the topology — frequency creates both DESIGN PRINCIPLES (from 4-checkpoint audit): CP1: Potential evolves WITH the pattern, not imposed beforehand CP2: Energy accumulation makes potential scalar (deepens with time) CP3: Chirality emerges from temporal progression (frame sequence) CP4: Frequency creates topology — interference pattern IS the surface MECHANISM: Each cycle: 1. Pulse: N=3 sources inject energy (f) 2. Rest: decoherence gap allows crystallization (t) 3. Accumulate: intensity adds to running potential V(x) 4. Feedback: V(x) modifies r(x) for the NEXT cycle 5. The system self-organizes: bright spots → deeper V → different r → different accumulation → structure emerges This is f+A|t: the amplitude (A) is the inverse of accumulated structure. More structure → lower A → changes the next pulse's interaction. MEASUREMENTS: Per-cycle snapshots of: - Accumulated intensity pattern I_acc(x,y) - Potential field V(x,y) = α * I_acc - Decoherence ratio field r(x,y) - Pattern CV (site differentiation) - Spatial frequency content (FFT) - N-fold symmetry analysis (angular autocorrelation) - Peak coordination structure - Chirality metric (rotation of intensity maxima between frames) DESIGNED FOR HETZNER: Progress after every snapshot, intermediate saves. CRITICAL DESIGN RULE: OUTCOME-AGNOSTIC. No pass/fail. The data shows what it shows. NO nudges: feedback coupling is a single parameter α, not tuned per-wave. Control run with α=0 (no feedback) provides baseline. """ import numpy as np from scipy import ndimage from scipy.spatial import cKDTree from collections import defaultdict import json import time import os def log(msg): print(msg, flush=True) # ====================================================================== # SECTION 1: SELF-CONSISTENT FDTD ENGINE # ====================================================================== def run_self_consistent_fdtd(grid_size=256, n_wavelengths=10, r_base=0.3, feedback_alpha=0.0, n_cycles=200, snapshot_interval=10, k=2*np.pi): """ 2D pulsed FDTD with self-consistent potential feedback. The wave equation: ∂²ψ/∂t² = c²∇²ψ + J(x,t; r(x)) where r(x) evolves each cycle based on accumulated intensity: V_n(x) = α * I_accumulated_n(x) (potential from energy coalescence) r_n(x) = r_base + V_n(x) (decoherence ratio from potential) r_n(x) clipped to [0.01, 0.49] Parameters: feedback_alpha: coupling strength. 0 = no feedback (control). Positive = bright regions get LONGER rest (more decoherence). This matches the theory: energy coalescence → curvature → increased decoherence space. n_cycles: total pulse-rest cycles snapshot_interval: save full state every N cycles """ L = n_wavelengths dx = L / grid_size c = 1.0 dt = 0.4 * dx / c # CFL with safety margin wavelength = 1.0 f0 = c / wavelength T = 1.0 / f0 steps_per_cycle = int(T / dt) # Grid x = np.linspace(-L/2, L/2, grid_size) y = np.linspace(-L/2, L/2, grid_size) X, Y = np.meshgrid(x, y) # N=3 source directions angles = np.array([0, 2*np.pi/3, 4*np.pi/3]) kx = k * np.cos(angles) ky = k * np.sin(angles) # Fields psi = np.zeros((grid_size, grid_size)) psi_prev = np.zeros((grid_size, grid_size)) # Accumulated intensity (the evolving potential) I_accumulated = np.zeros((grid_size, grid_size)) # Current decoherence ratio field r_field = np.full((grid_size, grid_size), r_base) courant2 = (c * dt / dx) ** 2 # Storage for snapshots snapshots = [] # Track peak positions for chirality analysis prev_peak_positions = None for cycle in range(n_cycles): # ---- PULSE-REST CYCLE ---- # During each cycle, use POSITION-DEPENDENT duty cycle. # Sites with higher r(x) get fewer active steps. # Implemented as a hard step function: source ON where # cycle_phase < (1 - r(x)), OFF otherwise. I_this_cycle = np.zeros((grid_size, grid_size)) for step in range(steps_per_cycle): t_local = step * dt cycle_phase = step / steps_per_cycle # Position-dependent source: ON where r(x) < threshold # Hard step: source ON where cycle_phase < (1 - r(x)), OFF otherwise source_mask = (cycle_phase < (1.0 - r_field)).astype(float) J = np.zeros((grid_size, grid_size)) for j in range(3): J += source_mask * np.sin(kx[j]*X + ky[j]*Y - 2*np.pi*f0*(cycle*T + t_local)) # FDTD update laplacian = np.zeros_like(psi) laplacian[1:-1, 1:-1] = ( psi[2:, 1:-1] + psi[:-2, 1:-1] + psi[1:-1, 2:] + psi[1:-1, :-2] - 4 * psi[1:-1, 1:-1] ) psi_next = 2*psi - psi_prev + courant2 * laplacian + dt**2 * J # Mur absorbing boundaries mur = (c*dt - dx) / (c*dt + dx) psi_next[:, 0] = psi[:, 1] + mur * (psi_next[:, 1] - psi[:, 0]) psi_next[:, -1] = psi[:, -2] + mur * (psi_next[:, -2] - psi[:, -1]) psi_next[0, :] = psi[1, :] + mur * (psi_next[1, :] - psi[0, :]) psi_next[-1, :] = psi[-2, :] + mur * (psi_next[-2, :] - psi[-1, :]) psi_prev = psi.copy() psi = psi_next.copy() # Accumulate this cycle's intensity I_this_cycle += psi**2 I_this_cycle /= steps_per_cycle # ---- FEEDBACK: Update accumulated intensity and potential ---- I_accumulated += I_this_cycle # Normalize accumulated intensity to prevent runaway I_norm = I_accumulated / (cycle + 1) # Update potential: V(x) = α * normalized accumulated intensity V_field = feedback_alpha * I_norm # Update decoherence ratio field r_field = np.clip(r_base + V_field, 0.01, 0.49) # ---- SNAPSHOT ---- if cycle % snapshot_interval == 0 or cycle == n_cycles - 1: snapshot = analyze_snapshot( I_norm, r_field, V_field, grid_size, cycle, prev_peak_positions ) snapshot["cycle"] = cycle snapshot["feedback_alpha"] = feedback_alpha # Update peak tracking for chirality if snapshot.get("peak_positions") is not None: prev_peak_positions = snapshot["peak_positions"] snapshots.append(snapshot) log(f" cycle {cycle:4d}: CV={snapshot['cv']:.5f} " f"nPeak={snapshot['n_peaks']:3d} " f"contrast={snapshot['contrast']:.3f} " f"r_range=[{snapshot['r_min']:.3f},{snapshot['r_max']:.3f}] " f"V_max={snapshot['V_max']:.4f} " f"5fold={snapshot.get('sym_5', 0):.3f} " f"chiral={snapshot.get('chirality_angle', 0):.2f}°") return snapshots, I_norm, r_field def analyze_snapshot(I_field, r_field, V_field, grid_size, cycle, prev_peaks=None): """Analyze current state for structure metrics.""" edge = int(grid_size * 0.15) interior = I_field[edge:-edge, edge:-edge] I_max = np.max(interior) I_min = np.min(interior) I_mean = np.mean(interior) if I_max < 1e-20: return { "cv": 0.0, "n_peaks": 0, "contrast": 0.0, "r_min": float(np.min(r_field)), "r_max": float(np.max(r_field)), "V_max": float(np.max(np.abs(V_field))), "sym_5": 0.0, "chirality_angle": 0.0, } contrast = (I_max - I_min) / (I_max + I_min) if (I_max + I_min) > 0 else 0 # Find peaks max_filt = ndimage.maximum_filter(interior, size=7) is_max = (interior == max_filt) threshold = I_min + 0.5 * (I_max - I_min) is_peak = is_max & (interior > threshold) peak_positions = np.argwhere(is_peak) peak_values = interior[is_peak] n_peaks = len(peak_values) cv = float(np.std(peak_values) / np.mean(peak_values)) if n_peaks > 1 else 0.0 # Angular symmetry analysis — check 3,4,5,6-fold for comparison sym_scores = {3: 0.0, 4: 0.0, 5: 0.0, 6: 0.0} if n_peaks >= 7: center = np.array([interior.shape[0]//2, interior.shape[1]//2]) relative = peak_positions - center angles_from_center = np.arctan2(relative[:, 0], relative[:, 1]) angle_diffs = np.abs(np.subtract.outer(angles_from_center, angles_from_center)) n_pairs = n_peaks * (n_peaks - 1) if n_pairs == 0: n_pairs = 1 # Random baseline: fraction of [0, π] within tolerance tolerance = 0.15 random_baseline = tolerance / np.pi for n_fold in [3, 4, 5, 6]: target = 2 * np.pi / n_fold # Check ALL harmonics of this fold (e.g., 72°, 144°, 216°, 288° for 5-fold) count = 0 for harmonic in range(1, n_fold): target_h = harmonic * 2 * np.pi / n_fold count += np.sum(np.abs(angle_diffs - target_h) < tolerance) # Normalize by number of pairs and harmonics, subtract random baseline raw = count / (n_pairs * (n_fold - 1)) sym_scores[n_fold] = float(max(0, raw - random_baseline)) sym_5 = sym_scores.get(5, 0.0) # Chirality: angular rotation of peak pattern between snapshots # Uses intensity-weighted centroid of ALL peaks (not top-N selection) # to avoid peak-set instability artifacts chirality_angle = 0.0 if prev_peaks is not None and len(prev_peaks) >= 3 and n_peaks >= 3: center = np.array([interior.shape[0]//2, interior.shape[1]//2]) # Current frame: intensity-weighted centroid of all peaks if len(peak_values) > 0: weights = peak_values / np.sum(peak_values) curr_centroid = np.average(peak_positions, axis=0, weights=weights) curr_angle = np.arctan2(curr_centroid[0]-center[0], curr_centroid[1]-center[1]) # Previous frame: uniform-weighted centroid (no intensity info stored) prev_arr = np.array(prev_peaks) prev_centroid = np.mean(prev_arr, axis=0) prev_angle = np.arctan2(prev_centroid[0]-center[0], prev_centroid[1]-center[1]) # Angle difference with wrapping to [-180, 180] diff = curr_angle - prev_angle diff = np.arctan2(np.sin(diff), np.cos(diff)) # proper wrapping chirality_angle = float(np.degrees(diff)) return { "cv": cv, "n_peaks": n_peaks, "contrast": float(contrast), "I_max": float(I_max), "I_mean": float(I_mean), "r_min": float(np.min(r_field)), "r_max": float(np.max(r_field)), "r_mean": float(np.mean(r_field)), "V_max": float(np.max(np.abs(V_field))), "sym_3": sym_scores.get(3, 0.0), "sym_4": sym_scores.get(4, 0.0), "sym_5": sym_5, "sym_6": sym_scores.get(6, 0.0), "chirality_angle": chirality_angle, "peak_positions": peak_positions.tolist() if n_peaks > 0 else None, } # ====================================================================== # MAIN # ====================================================================== def main(): log("=" * 70) log("B.6.7 SELF-CONSISTENT f+A|t UNFOLDING TEST") log("Does feedback-coupled pulsed interference self-organize?") log("=" * 70) log("") grid_size = 256 n_wavelengths = 10 n_cycles = 200 snapshot_interval = 10 r_base = 0.3 # Feedback strengths to test # α=0: control (no feedback, matches flat potential baseline) # α>0: energy coalescence deepens potential feedback_values = [0.0, 0.001, 0.005, 0.01, 0.02, 0.05, 0.1, 0.2] log(f"Grid: {grid_size}x{grid_size}, {n_wavelengths} wavelengths") log(f"Cycles: {n_cycles}, snapshots every {snapshot_interval}") log(f"r_base: {r_base}") log(f"Feedback α values: {feedback_values}") log("") all_results = {} t_start = time.time() for alpha in feedback_values: log("=" * 70) log(f"FEEDBACK α = {alpha}") if alpha == 0: log(" CONTROL: No feedback (flat potential, should match baseline)") else: log(f" Active: accumulated intensity → potential → decoherence") log("=" * 70) t0 = time.time() snapshots, I_final, r_final = run_self_consistent_fdtd( grid_size=grid_size, n_wavelengths=n_wavelengths, r_base=r_base, feedback_alpha=alpha, n_cycles=n_cycles, snapshot_interval=snapshot_interval, ) elapsed = time.time() - t0 log(f" Elapsed: {elapsed:.0f}s") log("") # Summary for this alpha if snapshots: first = snapshots[0] last = snapshots[-1] log(f" Evolution: CV {first['cv']:.5f} → {last['cv']:.5f} " f"r_range [{first['r_min']:.3f},{first['r_max']:.3f}] → " f"[{last['r_min']:.3f},{last['r_max']:.3f}] " f"V_max {first['V_max']:.4f} → {last['V_max']:.4f}") # Check if differentiation grew over time cv_growth = last['cv'] - first['cv'] if cv_growth > 0.01: log(f" ** DIFFERENTIATION GROWING: +{cv_growth:.5f} **") # Check if r_field spread increased (potential evolving) r_spread_first = first['r_max'] - first['r_min'] r_spread_last = last['r_max'] - last['r_min'] if r_spread_last > r_spread_first + 0.01: log(f" ** POTENTIAL DEEPENING: r spread {r_spread_first:.3f} → {r_spread_last:.3f} **") # Check for chirality (consistent rotation direction) chiral_angles = [s.get('chirality_angle', 0) for s in snapshots[1:]] if chiral_angles: mean_chiral = np.mean(chiral_angles) if abs(mean_chiral) > 1.0: direction = "RIGHT" if mean_chiral > 0 else "LEFT" log(f" ** CHIRALITY: {direction}-handed, mean rotation {mean_chiral:.2f}°/snapshot **") log("") # Strip peak_positions from snapshots for JSON (too large) clean_snapshots = [] for s in snapshots: cs = {k: v for k, v in s.items() if k != "peak_positions"} clean_snapshots.append(cs) all_results[f"alpha_{alpha}"] = { "feedback_alpha": alpha, "snapshots": clean_snapshots, "elapsed": elapsed, "r_base": r_base, } # Intermediate save _save_intermediate(all_results) total_elapsed = time.time() - t_start # ---- CROSS-ALPHA COMPARISON ---- log("=" * 70) log("CROSS-ALPHA COMPARISON") log("=" * 70) log("") log(f"{'α':>8s} | {'CV_0':>8s} {'CV_end':>8s} {'ΔCV':>8s} | " f"{'r_spread':>8s} {'V_max':>8s} | {'5fold':>6s} {'chiral':>7s}") log("-" * 75) for alpha in feedback_values: key = f"alpha_{alpha}" if key in all_results: snaps = all_results[key]["snapshots"] if snaps: first = snaps[0] last = snaps[-1] delta_cv = last['cv'] - first['cv'] r_spread = last['r_max'] - last['r_min'] chiral_angles = [s.get('chirality_angle', 0) for s in snaps[1:]] mean_chiral = np.mean(chiral_angles) if chiral_angles else 0 flag = "" if delta_cv > 0.01: flag = " ***" log(f"{alpha:8.3f} | {first['cv']:8.5f} {last['cv']:8.5f} " f"{delta_cv:8.5f} | {r_spread:8.4f} {last['V_max']:8.4f} | " f"{last.get('sym_5',0):6.3f} {mean_chiral:7.2f}{flag}") log("") log(f"Total elapsed: {total_elapsed:.0f}s") _save_final(all_results, total_elapsed) def _save_intermediate(results): output_dir = os.path.join( os.path.dirname(os.path.abspath(__file__)), "..", "..", "tlt notes", "theory", "mathematical_framework" ) path = os.path.join(output_dir, "B6_self_consistent_partial.json") try: with open(path, 'w') as f: json.dump(results, f, indent=2, default=_json_default) except Exception: pass def _save_final(results, elapsed): output_dir = os.path.join( os.path.dirname(os.path.abspath(__file__)), "..", "..", "tlt notes", "theory", "mathematical_framework" ) json_path = os.path.join(output_dir, "B6_self_consistent_results.json") txt_path = os.path.join(output_dir, "B6_self_consistent_results.txt") with open(json_path, 'w') as f: json.dump(results, f, indent=2, default=_json_default) with open(txt_path, 'w') as f: f.write("=" * 70 + "\n") f.write("B.6.7 SELF-CONSISTENT f+A|t UNFOLDING — RESULTS\n") f.write("=" * 70 + "\n") f.write(f"Date: {time.strftime('%Y-%m-%d %H:%M:%S')}\n") f.write(f"Elapsed: {elapsed:.0f}s\n\n") f.write("See JSON for full data.\n") log(f" Results: {txt_path}") log(f" JSON: {json_path}") def _json_default(obj): if isinstance(obj, (np.integer,)): return int(obj) if isinstance(obj, (np.floating,)): return float(obj) if isinstance(obj, np.ndarray): return obj.tolist() return str(obj) if __name__ == "__main__": main()