import gradio as gr
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# --- Constants ---
Z = 115  # Moscovium
A_BASE = 290  # Starting mass number

# SEMF constants (Weizsäcker)
a_v = 15.8
a_s = 18.3
a_c = 0.714
a_a = 23.2
a_p = 12.0


# --- Physics Functions ---

def pairing_term(A, Z):
    """Correct SEMF pairing term."""
    N = A - Z

    if A % 2 == 1:
        return 0.0  # odd-A

    elif Z % 2 == 0 and N % 2 == 0:
        return +a_p / np.sqrt(A)  # even-even

    else:
        return -a_p / np.sqrt(A)  # odd-odd


def binding_energy(A, Z):
    """Weizsäcker semi-empirical mass formula."""
    volume = a_v * A
    surface = a_s * A ** (2/3)
    coulomb = a_c * Z * (Z - 1) / A ** (1/3)
    asymmetry = a_a * (A - 2 * Z) ** 2 / A
    pairing = pairing_term(A, Z)

    return volume - surface - coulomb - asymmetry + pairing


def binding_energy_per_nucleon(E_bind, A):
    return E_bind / A


def stability_probability(BE_per_A):
    """
    Stability based on binding energy per nucleon.
    Peaks around ~7–8 MeV/A.
    """
    return 1 / (1 + np.exp(-2.0 * (BE_per_A - 7.5)))


def half_life_estimate(stability):
    """
    Stochastic half-life approximation.
    Adds noise to simulate chaotic decay behavior.
    """
    noise = np.random.normal(0, 0.2)
    return np.exp(4 * stability + noise)


def gravitational_effect(mass, speculative=False, coupling=1.0):
    """
    Mass-energy gravitational proxy.
    This mode introduces geometry-aware modulation.
    """
    G_effect = mass * 1e-27

    if speculative:
        density_factor = mass ** (1/3)
        resonance = np.sin(mass * 0.01)
        G_effect *= (1 + coupling * resonance / density_factor)

    return G_effect


# --- Simulation ---

def run_simulation(max_neutrons, speculative, coupling):
    results = []
    A = A_BASE

    for n in range(int(max_neutrons)):
        A += 1

        E_bind = binding_energy(A, Z)
        BE_per_A = binding_energy_per_nucleon(E_bind, A)
        stability = stability_probability(BE_per_A)
        half_life = half_life_estimate(stability)
        grav = gravitational_effect(A, speculative, coupling)

        results.append({
            "Neutrons Added": n + 1,
            "Mass Number": A,
            "Binding Energy (MeV)": E_bind,
            "BE/Nucleon (MeV)": BE_per_A,
            "Stability": stability,
            "Half-Life (arb)": half_life,
            "Gravitational Effect": grav,
        })

        # Improved break condition
        if BE_per_A < 6.5 and stability < 0.2:
            break

    df = pd.DataFrame(results)

    # --- Derived Signal (Second Derivative) ---
    if len(df) > 2:
        df["d2E"] = np.gradient(np.gradient(df["Binding Energy (MeV)"]))
    else:
        df["d2E"] = 0

    # --- Plotting ---
    fig, axs = plt.subplots(4, 1, figsize=(7, 13))
    fig.patch.set_facecolor("#0d0d0d")

    plot_cfg = [
        ("Binding Energy (MeV)", "Binding Energy"),
        ("BE/Nucleon (MeV)", "Binding Energy per Nucleon"),
        ("Stability", "Stability"),
        ("Gravitational Effect", "Gravitational Effect"),
    ]

    for ax, (col, title) in zip(axs, plot_cfg):
        ax.plot(df["Neutrons Added"], df[col], linewidth=1.8)
        ax.set_title(title, color="#e0e0e0", fontsize=11)
        ax.set_facecolor("#141414")
        ax.tick_params(colors="#888888")
        for spine in ax.spines.values():
            spine.set_edgecolor("#333333")

    plt.tight_layout(pad=2.0)

    return df.round(5), fig


# --- Gradio UI ---

with gr.Blocks(theme=gr.themes.Monochrome()) as demo:
    gr.Markdown("""
# 🧪 Vers3Dynamics Moscovium Neutron

This simulation explores nuclear stability trends for **Z = 115 (Moscovium)** using the
semi-empirical Weizsäcker mass formula.

A layer models hypothetical coupling between nuclear density and spacetime curvature.

**Stability boundary:** Binding energy per nucleon below ~6.5 MeV/A.

⚠️ The gravitational model is experimental.
""")

    with gr.Row():
        max_neutrons = gr.Slider(10, 120, value=88, step=1,
                                 label="Max Neutrons to Add")

        speculative = gr.Checkbox(label="Enable Gravitational Coupling")

        coupling = gr.Slider(0.1, 5.0, value=1.0, step=0.1,
                             label="Coupling Strength")

    run_btn = gr.Button("▶ Run Simulation", variant="primary")

    output_table = gr.Dataframe(label="Simulation Output")
    output_plot = gr.Plot(label="State Curves")

    run_btn.click(
        fn=run_simulation,
        inputs=[max_neutrons, speculative, coupling],
        outputs=[output_table, output_plot]
    )

demo.launch()