#!/usr/bin/env python3
"""Educational thermal-processing model summaries for SemiAgora replays.

The calculations below create public teaching datasets for anneal and thermal
budget intuition. They are normalized, precomputed replays. They are not tool
recipes, qualified process windows, wafer-release criteria, or measured
SemiAgora process data.
"""

from __future__ import annotations

import json
import math
from pathlib import Path

K_B_EV = 8.617333262145e-5
OUT = Path(__file__).resolve().parent


def source(label: str, url: str, note: str) -> dict[str, str]:
    return {"label": label, "url": url, "note": note}


SOURCES = {
    "mit6152_implant": source(
        "MIT OCW 6.152J ion implantation lecture",
        "https://ocw.mit.edu/courses/6-152j-micro-nano-processing-technology-fall-2005/resources/lecture6/",
        "Public lecture-note page anchoring diffusion, ion implantation, projected range, and implantation damage vocabulary.",
    ),
    "mit6774_course": source(
        "MIT OCW 6.774 front-end processing",
        "https://ocw.mit.edu/courses/6-774-physics-of-microfabrication-front-end-processing-fall-2004/",
        "Public course page covering oxidation, diffusion, ion implantation, epitaxy, and front-end process physics.",
    ),
    "mit6774_ted": source(
        "MIT OCW 6.774 transient enhanced diffusion",
        "https://ocw.mit.edu/courses/6-774-physics-of-microfabrication-front-end-processing-fall-2004/resources/mit6_774f04_lec15_mp4/",
        "Public lecture page for TED simulation examples and post-implant diffusion context.",
    ),
    "rta_control": source(
        "Optimal control of rapid thermal annealing",
        "https://web.mit.edu/braatzgroup/70_Optimal_control_of_rapid_thermal_annealing_in_a_semiconductor_process.pdf",
        "Open PDF discussing RTA, ultrashallow junctions, TED, sheet-resistance tradeoffs, and temperature trajectory control.",
    ),
    "comsol_rta": source(
        "COMSOL Rapid Thermal Annealing model",
        "https://www.comsol.com/model/rapid-thermal-annealing-504",
        "Public model page describing RTA as dopant activation and metal-contact interfacial reaction with transient heat transfer.",
    ),
    "mit6774_silicide": source(
        "MIT OCW 6.774 silicides and contacts",
        "https://ocw.mit.edu/courses/6-774-physics-of-microfabrication-front-end-processing-fall-2004/resources/mit6_774f04_lec22_mp4/",
        "Public lecture page anchoring silicide/contact vocabulary and front-end contact formation context.",
    ),
    "nist_facility": source(
        "NIST Semiconductor Electronics Division facility summary",
        "https://nvlpubs.nist.gov/nistpubs/Legacy/IR/nistir6430.pdf",
        "Public NIST report listing microfabrication facility capabilities including furnaces, thermal oxidation, CVD, and annealing.",
    ),
    "nist_laser": source(
        "NIST laser annealed arsenic ultrashallow junctions",
        "https://www.nist.gov/publications/deactivation-sub-melt-laser-annealed-arsenic-ultra-shallow-junctions-silicon-during",
        "Public NIST publication page for laser anneal, activation/deactivation, Hall, EXAFS, and SIMS context.",
    ),
}


COMMON_LIMITS = [
    "No furnace recipe, RTA program, laser scan condition, gas ambient, wafer handling, ramp-rate setpoint, endpoint, or tool-control instruction is provided.",
    "All axes are normalized teaching axes or simplified time-temperature profiles and must not be copied into equipment settings.",
    "The values are synthetic teaching replays, not measured SemiAgora wafer data, qualified metrology, or process release limits.",
    "Real anneal work requires material stack review, contamination controls, calibrated temperature metrology, electrical metrology, safety approval, and tool-specific training.",
]


def arrhenius_weight(temp_c: float, activation_ev: float, ref_c: float = 1000.0) -> float:
    temp_k = temp_c + 273.15
    ref_k = ref_c + 273.15
    return math.exp((-activation_ev / K_B_EV) * ((1.0 / temp_k) - (1.0 / ref_k)))


def weighted_seconds(profile: list[tuple[float, float]], activation_ev: float) -> float:
    return sum(seconds * arrhenius_weight(temp_c, activation_ev) for temp_c, seconds in profile)


def activation_fraction(act_budget_s: float, scale_s: float = 1.35, ceiling: float = 0.94) -> float:
    return ceiling * (1.0 - math.exp(-act_budget_s / scale_s))


def diffusion_length_nm(diff_budget_s: float, scale_nm: float = 6.25) -> float:
    return scale_nm * math.sqrt(max(diff_budget_s, 0.0))


def rounded(value: float, ndigits: int = 2) -> float:
    return round(value, ndigits)


def write_json(name: str, data: dict) -> None:
    (OUT / name).write_text(json.dumps(data, indent=2) + "\n", encoding="utf-8")


def build_activation_diffusion() -> dict:
    profiles = [
        {
            "id": "low_exposure",
            "label": "Low-exposure anneal",
            "profile_c_s": [(650, 6), (875, 5), (650, 6)],
            "teaching_note": "Preserves the shallow profile but leaves a large inactive fraction.",
        },
        {
            "id": "spike_rta",
            "label": "Spike RTA teaching profile",
            "profile_c_s": [(650, 3), (900, 2), (1050, 2), (850, 5)],
            "teaching_note": "Balances activation against diffusion broadening in the simplified replay.",
        },
        {
            "id": "long_soak",
            "label": "Longer soak teaching profile",
            "profile_c_s": [(800, 30), (1000, 30), (800, 30)],
            "teaching_note": "Improves activation and defect cleanup but broadens the junction proxy.",
        },
    ]
    cases = []
    for profile in profiles:
        act_budget = weighted_seconds(profile["profile_c_s"], 2.25)
        diff_budget = weighted_seconds(profile["profile_c_s"], 3.55)
        activation = activation_fraction(act_budget)
        diffusion = diffusion_length_nm(diff_budget)
        sheet_r = 285.0 / max(activation, 0.05) + diffusion * 4.6
        cases.append(
            {
                **profile,
                "effective_activation_budget_s_at_1000c": rounded(act_budget),
                "effective_diffusion_budget_s_at_1000c": rounded(diff_budget),
                "activated_fraction_percent": rounded(activation * 100.0),
                "diffusion_length_nm": rounded(diffusion),
                "junction_depth_proxy_nm": rounded(28.0 + diffusion * 1.2),
                "sheet_resistance_proxy_ohm_sq": rounded(sheet_r, 1),
                "inactive_fraction_percent": rounded((1.0 - activation) * 100.0),
            }
        )

    return {
        "schema": "semiagora.process-simulation.v1",
        "experiment_id": "SA-PROC-THERM-RTA-ACTIVATION-001",
        "title": "RTA activation and diffusion tradeoff",
        "execution_mode": "precomputed-only",
        "model_boundary": "Educational rapid-thermal-anneal replay using normalized Arrhenius budgets; not a recipe or measured wafer result.",
        "metric_contract": {
            "activation": "Higher activated_fraction_percent means more dopant is electrically useful in the teaching model.",
            "diffusion": "Higher diffusion_length_nm means more junction broadening risk in the teaching model.",
            "sheet_resistance": "Sheet resistance proxy combines activation benefit and diffusion penalty for comparison only.",
        },
        "cases": cases,
        "derived_metrics": {
            "best_balanced_case": "spike_rta",
            "low_exposure_activation_percent": cases[0]["activated_fraction_percent"],
            "spike_activation_percent": cases[1]["activated_fraction_percent"],
            "long_soak_junction_depth_proxy_nm": cases[2]["junction_depth_proxy_nm"],
        },
        "verification": {
            "case_count": len(cases),
            "activation_increases_with_budget": cases[0]["activated_fraction_percent"] < cases[1]["activated_fraction_percent"] < cases[2]["activated_fraction_percent"],
            "diffusion_increases_with_budget": cases[0]["diffusion_length_nm"] < cases[1]["diffusion_length_nm"] < cases[2]["diffusion_length_nm"],
            "all_cases_precomputed": True,
        },
        "limitations": COMMON_LIMITS,
        "sources": [SOURCES["mit6152_implant"], SOURCES["mit6774_course"], SOURCES["rta_control"], SOURCES["comsol_rta"]],
    }


def build_thermal_budget() -> dict:
    cases = [
        {
            "id": "furnace_soak",
            "label": "Furnace soak teaching mode",
            "relative_peak_temperature": "moderate",
            "relative_duration": "long",
            "activation_score": 0.88,
            "diffusion_broadening_nm": 28.5,
            "uniformity_score": 0.92,
            "temperature_control_risk": 0.18,
        },
        {
            "id": "spike_rta",
            "label": "Spike RTA teaching mode",
            "relative_peak_temperature": "high",
            "relative_duration": "short",
            "activation_score": 0.84,
            "diffusion_broadening_nm": 10.4,
            "uniformity_score": 0.78,
            "temperature_control_risk": 0.42,
        },
        {
            "id": "millisecond_anneal",
            "label": "Millisecond anneal teaching mode",
            "relative_peak_temperature": "very high local",
            "relative_duration": "very short",
            "activation_score": 0.74,
            "diffusion_broadening_nm": 3.2,
            "uniformity_score": 0.66,
            "temperature_control_risk": 0.71,
        },
    ]
    for case in cases:
        case["thermal_budget_index"] = rounded(case["diffusion_broadening_nm"] / 28.5)
        case["shallow_junction_score"] = rounded(1.0 / (1.0 + case["diffusion_broadening_nm"] / 12.0))
        case["paper_review_warning"] = (
            "Ask for temperature metrology, spatial uniformity, and electrical activation evidence before comparing anneal modes."
        )

    return {
        "schema": "semiagora.process-simulation.v1",
        "experiment_id": "SA-PROC-THERM-BUDGET-COMPARE-001",
        "title": "Spike versus furnace thermal budget",
        "execution_mode": "precomputed-only",
        "model_boundary": "Educational comparison of normalized anneal modes; not an equipment program or process recommendation.",
        "cases": cases,
        "derived_metrics": {
            "lowest_diffusion_case": "millisecond_anneal",
            "highest_uniformity_case": "furnace_soak",
            "balanced_public_teaching_case": "spike_rta",
        },
        "verification": {
            "case_count": len(cases),
            "diffusion_decreases_with_shorter_duration": cases[0]["diffusion_broadening_nm"] > cases[1]["diffusion_broadening_nm"] > cases[2]["diffusion_broadening_nm"],
            "temperature_control_risk_increases_with_aggressive_mode": cases[0]["temperature_control_risk"] < cases[1]["temperature_control_risk"] < cases[2]["temperature_control_risk"],
            "all_cases_precomputed": True,
        },
        "limitations": COMMON_LIMITS,
        "sources": [SOURCES["rta_control"], SOURCES["comsol_rta"], SOURCES["nist_laser"], SOURCES["nist_facility"]],
    }


def build_damage_recovery() -> dict:
    times = [0, 1, 2, 5, 10, 20, 40]
    samples = []
    for t in times:
        defect_remaining = 100.0 * math.exp(-t / 8.5)
        activated = 90.0 * (1.0 - math.exp(-t / 6.0))
        ted_tail = 5.0 + 28.0 * (1.0 - math.exp(-t / 18.0))
        samples.append(
            {
                "normalized_anneal_time": t,
                "defect_remaining_percent": rounded(defect_remaining),
                "activated_fraction_percent": rounded(activated),
                "ted_tail_proxy_nm": rounded(ted_tail),
                "review_note": "Damage repair, activation, and TED tail growth must be reviewed together.",
            }
        )

    return {
        "schema": "semiagora.process-simulation.v1",
        "experiment_id": "SA-PROC-THERM-DAMAGE-RECOVERY-001",
        "title": "Implant damage recovery and TED",
        "execution_mode": "precomputed-only",
        "model_boundary": "Educational post-implant recovery replay with normalized anneal time; not a calibrated defect or diffusion model.",
        "samples": samples,
        "derived_metrics": {
            "half_damage_recovery_time_index": 6.0,
            "activation_at_time_10_percent": next(x["activated_fraction_percent"] for x in samples if x["normalized_anneal_time"] == 10),
            "ted_tail_at_time_40_nm": samples[-1]["ted_tail_proxy_nm"],
        },
        "verification": {
            "sample_count": len(samples),
            "defects_decrease": all(samples[i]["defect_remaining_percent"] >= samples[i + 1]["defect_remaining_percent"] for i in range(len(samples) - 1)),
            "activation_increases": all(samples[i]["activated_fraction_percent"] <= samples[i + 1]["activated_fraction_percent"] for i in range(len(samples) - 1)),
            "ted_tail_increases": all(samples[i]["ted_tail_proxy_nm"] <= samples[i + 1]["ted_tail_proxy_nm"] for i in range(len(samples) - 1)),
            "all_cases_precomputed": True,
        },
        "limitations": COMMON_LIMITS + [
            "Point-defect clustering, dopant species, channeling, amorphization, surface recombination, and material-specific calibration are not modeled."
        ],
        "sources": [SOURCES["mit6152_implant"], SOURCES["mit6774_ted"], SOURCES["rta_control"]],
    }


def build_silicide_window() -> dict:
    cases = [
        {
            "id": "under_reacted",
            "label": "Under-reacted contact",
            "normalized_anneal_index": 0.25,
            "silicide_completion_percent": 38,
            "sheet_resistance_ohm_sq": 18.5,
            "agglomeration_risk": 0.04,
            "phase_note": "Incomplete reaction dominates resistance.",
        },
        {
            "id": "target_window",
            "label": "Target teaching window",
            "normalized_anneal_index": 0.58,
            "silicide_completion_percent": 88,
            "sheet_resistance_ohm_sq": 4.4,
            "agglomeration_risk": 0.16,
            "phase_note": "Low-resistance contact region in the teaching replay.",
        },
        {
            "id": "over_annealed",
            "label": "Over-annealed contact",
            "normalized_anneal_index": 0.9,
            "silicide_completion_percent": 95,
            "sheet_resistance_ohm_sq": 9.8,
            "agglomeration_risk": 0.72,
            "phase_note": "Agglomeration or morphology risk dominates the non-claim boundary.",
        },
    ]
    return {
        "schema": "semiagora.process-simulation.v1",
        "experiment_id": "SA-PROC-THERM-SILICIDE-CONTACT-001",
        "title": "Silicide contact anneal window",
        "execution_mode": "precomputed-only",
        "model_boundary": "Educational silicide/contact anneal window replay; not a metal stack recipe, phase diagram, or contact qualification.",
        "cases": cases,
        "derived_metrics": {
            "lowest_resistance_case": "target_window",
            "target_sheet_resistance_ohm_sq": 4.4,
            "over_anneal_agglomeration_risk": 0.72,
        },
        "verification": {
            "case_count": len(cases),
            "minimum_resistance_case_is_target": min(cases, key=lambda x: x["sheet_resistance_ohm_sq"])["id"] == "target_window",
            "agglomeration_risk_increases": cases[0]["agglomeration_risk"] < cases[1]["agglomeration_risk"] < cases[2]["agglomeration_risk"],
            "all_cases_precomputed": True,
        },
        "limitations": COMMON_LIMITS + [
            "Specific metals, phase transitions, line width effects, silicon consumption, stress, capping layers, and contact-chain metrology are not modeled."
        ],
        "sources": [SOURCES["mit6774_silicide"], SOURCES["comsol_rta"], SOURCES["nist_facility"]],
    }


def main() -> None:
    datasets = {
        "process-rta-activation-diffusion-web-v1.json": build_activation_diffusion(),
        "process-spike-vs-furnace-thermal-budget-web-v1.json": build_thermal_budget(),
        "process-implant-damage-recovery-web-v1.json": build_damage_recovery(),
        "process-silicide-contact-anneal-window-web-v1.json": build_silicide_window(),
    }
    for name, data in datasets.items():
        write_json(name, data)
    print(json.dumps({name: data["experiment_id"] for name, data in datasets.items()}, indent=2))


if __name__ == "__main__":
    main()
