import numpy as np

phi = (1 + np.sqrt(5)) / 2
PRIMES = [
    2, 3, 5, 7, 11, 13, 17, 19, 23, 29,
    31, 37, 41, 43, 47, 53, 59, 61, 67, 71,
    73, 79, 83, 89, 97, 101, 103, 107, 109, 113,
    127, 131, 137, 139, 149, 151, 157, 163, 167, 173,
    179, 181, 191, 193, 197, 199, 211, 223, 227, 229
]

def fib_real(n):
    from math import cos, pi, sqrt
    phi_inv = 1 / phi
    return (phi**n / sqrt(5)) - (phi_inv**n * cos(pi * n))

def Omega(n, alpha=0.1):
    return np.exp(-alpha * n)

def D_n(n, r=1.0, k=1.0, alpha=0.1):
    Fn = fib_real(n)
    idx = int(np.floor(n)) % len(PRIMES)
    Pn = PRIMES[idx]
    val = phi * Fn * (2 ** n) * Pn * Omega(n, alpha)
    val = max(val, 1e-15)
    return np.sqrt(val) * (r ** k)

# Manually select reference recursion depth where model speed ~ 1
n0 = 4.0  # example, corresponds to your previous intuition
c_ref = D_n(n0)

# Calculate scale factor to physical speed of light in km/s at n0
c_phys_actual = 299792.458
S = c_phys_actual / c_ref

def physical_speed_of_light(n):
    c_model = D_n(n)
    return c_model * S

import matplotlib.pyplot as plt

n_vals = np.linspace(0, 10, 100)
c_vals = [physical_speed_of_light(n) for n in n_vals]

plt.figure(figsize=(8,5))
plt.plot(n_vals, c_vals / 1e3)  # speed in thousands km/s for readability
plt.axhline(299.792, color='r', linestyle='--', label='Local c (km/s)')
plt.xlabel('Recursion Depth n')
plt.ylabel('Speed of Light (×10³ km/s)')
plt.title('Emergent Speed of Light vs. Recursion Depth')
plt.legend()
plt.grid(True)
plt.show()


print(f"Reference n0 = {n0}, D_n0 = {c_ref:.6f}, scale factor S = {S:.2f} km/s")

print("\nSpeed of light at various recursion depths:")
for n_test in np.linspace(0, 10, 11):
    c_phys = physical_speed_of_light(n_test)
    print(f"n = {n_test:4.2f}, c = {c_phys:12.2f} km/s")
