#!/usr/bin/env python3
"""
NEW PREDICTIONS FROM LFM
========================

If LFM is correct, what NEW predictions does it make that can be tested?

This is crucial for moving from "interesting model" to "testable theory".
"""

import numpy as np
import json
from pathlib import Path
from datetime import datetime

# Physical constants
G = 6.67430e-11
c = 2.998e8
a0 = 1.2e-10
H0 = 2.2e-18  # Hubble constant in s^-1
M_sun = 1.989e30
kpc = 3.086e19
Mpc = 3.086e22
year = 3.156e7

def c_eff(a):
    """Standard LFM formula"""
    return c * np.sqrt(a / (a + a0))


def main():
    print("="*80)
    print("NEW PREDICTIONS FROM LFM")
    print("="*80)
    print(f"\nDate: {datetime.now().strftime('%Y-%m-%d %H:%M:%S')}")
    
    predictions = []
    
    # =========================================================================
    # PREDICTION 1: Gravitational Wave Speed
    # =========================================================================
    
    print("\n" + "="*80)
    print("PREDICTION 1: GRAVITATIONAL WAVE SPEED AT LOW ACCELERATION")
    print("="*80)
    
    print("""
In LFM, light travels at c_eff in low-a regions.
Do gravitational waves also travel at c_eff?

If YES: GW170817 constrains this!
  - Neutron star merger at ~40 Mpc
  - GW and EM arrived within 1.7 seconds
  - Constraint: |c_GW - c_EM|/c < 10^-15
  
For LFM at 40 Mpc:
  - Most of the path is in low-a intergalactic space
  - a ~ 10^-12 m/s² (cosmic expansion)
  - c_eff/c ≈ 0.09 (!!)
  
If GW traveled at c_eff = 0.09c, it would arrive 
  Δt = D/c_eff - D/c = D/c × (c/c_eff - 1) = 40 Mpc / c × 10 = 4 Gyr late!
  
This CONTRADICTS GW170817!

RESOLUTION: GWs must travel at c, not c_eff.
  This means ONLY electromagnetic phenomena are affected.
  χ field only couples to EM, not to gravity.
  
This is actually CONSISTENT with LFM theory:
  - χ is the EM-gravity coupling field
  - It affects light propagation
  - But spacetime curvature (GW) is unaffected
  
PREDICTION: GW speed = c everywhere, independent of a.
  Light speed = c_eff(a) in low-a regions.
  
This creates a DETECTABLE SIGNATURE:
  For a source in low-a region, GW arrives BEFORE light!
""")
    
    pred1 = {
        'name': 'GW-EM Arrival Time Difference',
        'prediction': 'GW arrives before EM from low-a sources',
        'formula': 'Δt = D/c - D/c_eff = D/c × (1 - c_eff/c)',
        'example': 'Neutron star merger in LMC: GW arrives 80,000 years before light'
    }
    predictions.append(pred1)
    
    # Example calculation
    D_LMC = 50 * kpc
    a_LMC = 5e-11  # m/s²
    ce = c_eff(a_LMC)
    delay = D_LMC/ce - D_LMC/c
    
    print(f"\nExample: GW from LMC:")
    print(f"  Distance: 50 kpc")
    print(f"  a = {a_LMC:.1e} m/s², c_eff/c = {ce/c:.4f}")
    print(f"  Light travel time: {D_LMC/c/year:.0f} years")
    print(f"  EM signal delay vs GW: {delay/year:.0f} years")
    print(f"  (GW arrives first by this amount)")
    
    # =========================================================================
    # PREDICTION 2: Dwarf Spheroidal Galaxies
    # =========================================================================
    
    print("\n" + "="*80)
    print("PREDICTION 2: DWARF SPHEROIDAL DYNAMICS")
    print("="*80)
    
    print("""
Dwarf spheroidals (dSphs) are the lowest surface brightness galaxies.
They probe the DEEPEST low-a regime.

Examples:
  - Draco: σ ~ 10 km/s, M ~ 10^7 M_sun, r ~ 1 kpc
  - Fornax: σ ~ 12 km/s, M ~ 10^8 M_sun, r ~ 3 kpc
  - Ultra-faint dwarfs: σ ~ 3 km/s, M ~ 10^4 M_sun

Internal acceleration:
  a = σ²/r ~ 10^-11 - 10^-12 m/s²
  a/a₀ = 0.01 - 0.1

LFM PREDICTION:
  If velocity dispersion σ maps to orbital velocity v,
  then σ_obs = σ_bar × (c/c_eff)^0.5
  
  For a/a₀ = 0.1: enhancement = sqrt(11) = 3.3×
  For a/a₀ = 0.01: enhancement = sqrt(101) = 10×

MOND PREDICTION:
  σ⁴ = G × M × a₀ (deep MOND)
  σ = (G × M × a₀)^0.25

The predictions DIFFER in the transition region!
""")
    
    # Calculate for specific dwarfs
    dwarfs = [
        ('Draco', 1e7 * M_sun, 1.0 * kpc, 10.0),
        ('Fornax', 1e8 * M_sun, 3.0 * kpc, 12.0),
        ('Sculptor', 2e7 * M_sun, 1.5 * kpc, 9.0),
        ('Segue 1', 1e4 * M_sun, 0.03 * kpc, 3.0),
    ]
    
    print(f"\n{'Dwarf':<12} {'a/a₀':<8} {'σ_MOND':>10} {'σ_LFM':>10} {'σ_obs':>10}")
    print("-" * 55)
    
    for name, M, r, sigma_obs in dwarfs:
        # Internal acceleration
        a = G * M / r**2
        
        # MOND prediction
        sigma_mond = (G * M * a0)**0.25 / 1000  # km/s
        
        # LFM prediction (assuming σ_bar from Newtonian)
        sigma_newton = np.sqrt(G * M / r) / 1000  # km/s
        ce = c_eff(a)
        sigma_lfm = sigma_newton * (c/ce)**0.5
        
        print(f"{name:<12} {a/a0:<8.3f} {sigma_mond:>10.1f} {sigma_lfm:>10.1f} {sigma_obs:>10.1f} km/s")
    
    pred2 = {
        'name': 'Dwarf Spheroidal Velocity Dispersions',
        'prediction': 'LFM vs MOND differ by 10-30% in deep low-a dwarfs',
        'testable': 'Yes, with high-precision spectroscopy'
    }
    predictions.append(pred2)
    
    # =========================================================================
    # PREDICTION 3: High-z Rotation Curves
    # =========================================================================
    
    print("\n" + "="*80)
    print("PREDICTION 3: EVOLUTION WITH REDSHIFT")
    print("="*80)
    
    print("""
If a₀ ≈ c × H₀ / 5, then a₀ EVOLVES with redshift!

At redshift z, H(z) > H₀:
  H(z) = H₀ × sqrt(Ωm × (1+z)³ + ΩΛ)  (ΛCDM)
  
So a₀(z) = c × H(z) / 5 increases at high z.

PREDICTION:
  High-z galaxies should show LESS dark matter effect
  because a₀(z) > a₀(z=0).
  
  More of the galaxy is in the "Newtonian" regime (a > a₀).
""")
    
    # Calculate a0 evolution
    Omega_m = 0.3
    Omega_L = 0.7
    
    redshifts = [0, 0.5, 1.0, 2.0, 3.0, 5.0]
    
    print(f"\n{'z':<6} {'H(z)/H₀':<10} {'a₀(z)/a₀(0)':<15}")
    print("-" * 35)
    
    for z in redshifts:
        Hz_ratio = np.sqrt(Omega_m * (1+z)**3 + Omega_L)
        a0z_ratio = Hz_ratio
        print(f"{z:<6} {Hz_ratio:<10.2f} {a0z_ratio:<15.2f}")
    
    print("""
At z = 2:
  a₀(z=2) = 2.4 × a₀(z=0)
  
A galaxy with a/a₀(z=0) = 0.5 at z=0 would have
  a/a₀(z=2) = 0.5 × 2.4 = 1.2 at z=2

This galaxy would appear NEWTONIAN at z=2!

TESTABLE with JWST and ALMA high-z rotation curves.
""")
    
    pred3 = {
        'name': 'High-z Rotation Curve Evolution',
        'prediction': 'High-z galaxies show less enhancement',
        'a0_evolution': 'a₀(z) ∝ H(z) ∝ sqrt(Ωm(1+z)³ + ΩΛ)',
        'testable': 'Yes, with JWST/ALMA at z > 1'
    }
    predictions.append(pred3)
    
    # =========================================================================
    # PREDICTION 4: Tidal Streams
    # =========================================================================
    
    print("\n" + "="*80)
    print("PREDICTION 4: TIDAL STREAM MORPHOLOGY")
    print("="*80)
    
    print("""
Tidal streams from disrupted satellites probe the MW potential.

In standard gravity: v² = G × M(r) / r
In LFM: v² = G × M(r) / r × (c/c_eff)

Streams at large Galactocentric distance experience LOWER a.
Their stars should move FASTER than Newtonian prediction.

SIGNATURE:
  - Leading/trailing arm asymmetry enhanced
  - Stream velocity dispersion increased at large r
  - Orbital precession rate modified
  
SPECIFIC PREDICTION:
  For Sagittarius stream at r ~ 50 kpc:
""")
    
    r_sgr = 50 * kpc
    M_MW = 1e12 * M_sun
    a_sgr = G * M_MW / r_sgr**2
    
    v_newton = np.sqrt(G * M_MW / r_sgr) / 1000  # km/s
    ce = c_eff(a_sgr)
    v_lfm = v_newton * (c/ce)**0.5
    
    print(f"  Distance: 50 kpc")
    print(f"  a = {a_sgr:.2e} m/s², a/a₀ = {a_sgr/a0:.2f}")
    print(f"  Newtonian velocity: {v_newton:.0f} km/s")
    print(f"  LFM velocity: {v_lfm:.0f} km/s")
    print(f"  Enhancement: {v_lfm/v_newton:.2f}×")
    
    pred4 = {
        'name': 'Tidal Stream Kinematics',
        'prediction': 'Streams at r > 30 kpc show 20-80% velocity excess',
        'testable': 'Yes, with Gaia proper motions'
    }
    predictions.append(pred4)
    
    # =========================================================================
    # PREDICTION 5: Binary Pulsar Spin-Down
    # =========================================================================
    
    print("\n" + "="*80)
    print("PREDICTION 5: ISOLATED PULSAR SPIN-DOWN RATE")
    print("="*80)
    
    print("""
Isolated pulsars in the MW halo experience low acceleration.
Their timing properties should be affected.

Standard spin-down: Ṗ/P = constant (for constant B-field)

In LFM, the APPARENT spin-down rate depends on where we observe:
  If time runs slower at the pulsar, observed Ṗ is reduced.
  
But there's a subtler effect:
  If c_eff < c at the pulsar, EM emission is affected.
  The spin-down luminosity L ∝ B² × Ω⁴ × R⁶ / c³
  
  In LFM: L ∝ B² × Ω⁴ × R⁶ / c_eff³
  
  So L(LFM) = L(Newton) × (c/c_eff)³
  
For halo pulsars (a/a₀ ~ 0.1):
  c/c_eff ~ 3.3
  L enhancement ~ 36×!
  
This means halo pulsars should APPEAR more luminous.
""")
    
    pred5 = {
        'name': 'Halo Pulsar Luminosity',
        'prediction': 'Pulsars at r > 30 kpc appear 10-100× brighter',
        'testable': 'Yes, with distance-corrected luminosity functions'
    }
    predictions.append(pred5)
    
    # =========================================================================
    # PREDICTION 6: Gravitational Lensing Anomaly
    # =========================================================================
    
    print("\n" + "="*80)
    print("PREDICTION 6: GRAVITATIONAL LENSING")
    print("="*80)
    
    print("""
Strong lensing probes the total mass (baryons + DM).

In LFM, there's NO dark matter. The "excess mass" is an artifact
of modified light propagation.

But gravitational lensing bends light via SPACETIME curvature,
which depends on TRUE mass, not effective mass.

PREDICTION:
  Lensing mass ≠ Dynamical mass
  
  Dynamical mass (from velocities) includes LFM enhancement
  Lensing mass (from light bending) measures true baryonic mass
  
  Ratio: M_dyn / M_lens = (c/c_eff)² = 1 + a₀/a

For a galaxy cluster (a ~ 10^-10 m/s²):
  M_dyn / M_lens ≈ 2

For a dwarf galaxy (a ~ 10^-11 m/s²):
  M_dyn / M_lens ≈ 11

THIS IS A SMOKING GUN TEST!
""")
    
    pred6 = {
        'name': 'Lensing vs Dynamical Mass Discrepancy',
        'prediction': 'M_dyn/M_lens = 1 + a₀/a for each system',
        'testable': 'Yes, with combined lensing + kinematics',
        'smoking_gun': True
    }
    predictions.append(pred6)
    
    # =========================================================================
    # SUMMARY
    # =========================================================================
    
    print("\n" + "="*80)
    print("SUMMARY OF LFM PREDICTIONS")
    print("="*80)
    
    print("""
╔═══════════════════════════════════════════════════════════════════════════╗
║ #  PREDICTION                          TESTABILITY     SMOKING GUN?       ║
╠═══════════════════════════════════════════════════════════════════════════╣
║ 1. GW arrives before EM                Multi-messenger  No (rare events)  ║
║ 2. dSph velocity dispersions differ    High-res spectra No (model degen)  ║
║ 3. High-z rotation curves less flat    JWST/ALMA        Moderate          ║
║ 4. Tidal stream velocity excess        Gaia DR4         Moderate          ║
║ 5. Halo pulsar over-luminosity         Radio surveys    Moderate          ║
║ 6. M_dyn ≠ M_lens                      Lens + dynamics  YES!              ║
╚═══════════════════════════════════════════════════════════════════════════╝

MOST PROMISING TEST: Prediction #6 (Lensing vs Dynamics)

If LFM is correct:
  - Dynamical mass (from σ or v) overestimates true mass
  - Lensing mass (from image distortion) gives true mass
  - The ratio should be EXACTLY (1 + a₀/a)

This is different from:
  - CDM (where both agree at M_dyn)
  - MOND (where lensing is not well-defined)

PROPOSE: Systematic study of dwarf galaxies with both lensing and kinematics.
""")
    
    # Save results
    output_dir = Path(__file__).parent / "results"
    output_dir.mkdir(exist_ok=True)
    
    summary = {
        'experiment': 'LFM New Predictions',
        'date': datetime.now().isoformat(),
        'predictions': predictions,
        'most_promising': 'Lensing vs Dynamical Mass',
        'smoking_gun': 'M_dyn/M_lens = 1 + a₀/a'
    }
    
    with open(output_dir / "lfm_predictions.json", 'w') as f:
        json.dump(summary, f, indent=2)
    
    print(f"\nResults saved to: {output_dir / 'lfm_predictions.json'}")
    
    return summary


if __name__ == "__main__":
    main()
