Reproducing Wang, Hu and Irene 1994: thermal oxidation of germanium at 400 to 550 °C

One of the [studies we reproduced](/writing/interface-studies-reproduced) for [Questions we need Emily to answer](/writing/questions-for-emily-testing-a-theory-of-the-interface), the companion to [the essay](/writing/testing-a-theory-of-the-interface-with-an-ontology). Generated by `scripts/reproduce_studies.py` from `reproductions/wang_1994.py`. Every value from the paper was read by an AI agent and has not been checked by a person; the sample and history were re-read from the PDF alone by a second agent that had not seen our earlier transcription.

The study and why we depend on it

Y. Wang, Y. Z. Hu and E. A. Irene, “Electron cyclotron resonance plasma and thermal oxidation mechanisms of germanium”, ONR Technical Report 63 (1994), DTIC 94-09193; submitted to J. Vac. Sci. Technol.

The only germanium oxidation kinetics we hold, and so the only basis for the challenge's Ge/GeO2 system. A structure search over growth laws (scripts/search_growth.py) chose a power law in thickness on the 400, 450 and 550 °C points of Figure 1(b) and holds 500 °C out. This page writes the chamber sequence the paper describes around those points.

How well can this study be reproduced?

Moderate. The chamber sequence is stated: heat in vacuum, valve off, start ellipsometry, admit oxygen to 1 atm, evacuate, cool. What a simulation must supply: the heat-up time, how long the chamber takes to fill, the starting oxide (0 nm in Figure 1a, 1 nm in Figure 1b), and whether the long-time points came from this chamber or from a “conventional oxidation furnace”. Thicknesses are on a scanned figure; three are confirmed by a printed table, and a second printed table disagrees with the figure.

The grade is ours and is about the paper as a source for simulation, not about the quality of the experiment: high means the sample, every step and its timing are stated well enough that a simulation needs no choices of ours; moderate that some choices are ours and are listed; low that the result depends on choices the paper does not constrain.

The sample and its history, played back

The sample is drawn from the simulation: layer heights and the read-outs are what the ISL programs compute at each moment of the history. Steps the language cannot state (a clean, a metallization, a measurement) are on the clock and marked as recorded, not simulated.

Does our sample and history match the paper's?

Left: the paper's own words, with the page of the PDF. Right: what we simulate. 6 same, 4 approximated, 1 assumed, 4 not in ours.

The timing of the history

Time is counted from the start of the first simulated step. A step the language cannot state takes no simulated time; where the paper gives its duration, the duration is in the last column. Stated means the paper gives the duration; partly stated that it gives a bound or a range; not stated; assumed that the duration in our history is ours.

(100) Ge, 400 °C, 1 atm O2

Simulated time in this history: 9 h 51 min 30 s.

Across the samples

Timing the paper does not give

• Heat-up rate and any hold in vacuum before oxygen.
• How long the valved-off chamber took to reach 1 atm at 150 sccm (its volume is not given), while the clock was already running.
• Whether each point of Figure 1b is a separate sample or one sample cooled, measured and re-oxidized.
• Whether the long-time points came from this chamber or from a conventional furnace, and if so, how that furnace was loaded and ramped.
• Cooling rate and time in vacuum; time in air before any ex-situ measurement (XPS, FTIR, C-V).

What the second reading changed

• The stored recipe for each temperature is one hold in dry oxygen for the longest time. The paper's sequence has a heat-up in vacuum, a fill, the oxidation, an evacuation and a cooling in vacuum; written out, the heat-up and the cooling add about 0.04 nm in the program because it does not read the gas.
• The starting oxide is a fitted parameter of the program (“starting oxide fixed” in the selected structure, a value fitted with the others). The paper's two figures disagree on it (about 0 nm in Figure 1a, about 1 nm in Figure 1b); the sample cannot carry it.
• Table 3 prints 3 nm for an XPS sample at 400 °C, 540 min, where the dataset (Figure 1b) holds 4.4 nm at 556 min; and 8 nm at 550 °C, 10 min, where Figure 1a shows about 10 nm. Neither is in the dataset; both are compared here.
• The paper's time zero is the admission of oxygen into a closed chamber filling at 150 sccm; our curves start after a 30 s fill step, an offset of up to 5 min in the range searched.
• The program the predictor runs has no file: it is generated by scripts/search_growth.py from the selected structure, so this page cannot print it as ISL.

Their terms and ours

The simulation in our language

Each sample is one history file: a SAMPLE and a CONTROL_TRAJECTORY. This is the whole input; run it with PYTHONPATH=. python -m isl_harness.predict <file>.

(100) Ge, 400 °C, 1 atm O2

# (100) Ge, 400 °C, 1 atm O2
SAMPLE Ge100_400 {
    FAMILY: Ge_GeO2
    ORIENTATION: Ge100
    DOPING: n_type
    RESISTIVITY_OHM_CM: 20
    CLEAN: not_stated
}
CONTROL_TRAJECTORY History {
    EVENT Clean KIND clean TIMING assumed                        # Bromine-methanol etch, buffered HF, water rinse
    EVENT Load KIND load TIMING assumed                          # Loaded into the vacuum chamber
    RUN HeatUp DURATION 300 TEMP 25->400 AMBIENT vacuum TIMING bound VARY duration 60..1800
    EVENT Valve KIND measurement TIMING assumed                  # Pump valved off; real-time ellipsometry started
    RUN Fill DURATION 30 TEMP 400->400 AMBIENT mix:O2=0.5 TIMING bound VARY duration 5..300
    RUN Oxidation DURATION 33360 TEMP 400->400 AMBIENT dry_O2 TIMING bound
    EVENT Evacuate KIND other TIMING assumed                     # Chamber evacuated
    RUN Cool DURATION 1800 TEMP 400->25 AMBIENT vacuum PROFILE exponential TAU 300 TIMING bound
    EVENT SE KIND measurement TIMING assumed                     # Spectroscopic ellipsometry at room temperature, two-layer model
}

(100) Ge, 450 °C, 1 atm O2

# (100) Ge, 450 °C, 1 atm O2
SAMPLE Ge100_450 {
    FAMILY: Ge_GeO2
    ORIENTATION: Ge100
    DOPING: n_type
    RESISTIVITY_OHM_CM: 20
    CLEAN: not_stated
}
CONTROL_TRAJECTORY History {
    EVENT Clean KIND clean TIMING assumed                        # Bromine-methanol etch, buffered HF, water rinse
    EVENT Load KIND load TIMING assumed                          # Loaded into the vacuum chamber
    RUN HeatUp DURATION 300 TEMP 25->450 AMBIENT vacuum TIMING bound VARY duration 60..1800
    EVENT Valve KIND measurement TIMING assumed                  # Pump valved off; real-time ellipsometry started
    RUN Fill DURATION 30 TEMP 450->450 AMBIENT mix:O2=0.5 TIMING bound VARY duration 5..300
    RUN Oxidation DURATION 35820 TEMP 450->450 AMBIENT dry_O2 TIMING bound
    EVENT Evacuate KIND other TIMING assumed                     # Chamber evacuated
    RUN Cool DURATION 1800 TEMP 450->25 AMBIENT vacuum PROFILE exponential TAU 300 TIMING bound
    EVENT SE KIND measurement TIMING assumed                     # Spectroscopic ellipsometry at room temperature, two-layer model
}

(100) Ge, 500 °C, 1 atm O2

# (100) Ge, 500 °C, 1 atm O2
SAMPLE Ge100_500 {
    FAMILY: Ge_GeO2
    ORIENTATION: Ge100
    DOPING: n_type
    RESISTIVITY_OHM_CM: 20
    CLEAN: not_stated
}
CONTROL_TRAJECTORY History {
    EVENT Clean KIND clean TIMING assumed                        # Bromine-methanol etch, buffered HF, water rinse
    EVENT Load KIND load TIMING assumed                          # Loaded into the vacuum chamber
    RUN HeatUp DURATION 300 TEMP 25->500 AMBIENT vacuum TIMING bound VARY duration 60..1800
    EVENT Valve KIND measurement TIMING assumed                  # Pump valved off; real-time ellipsometry started
    RUN Fill DURATION 30 TEMP 500->500 AMBIENT mix:O2=0.5 TIMING bound VARY duration 5..300
    RUN Oxidation DURATION 35820 TEMP 500->500 AMBIENT dry_O2 TIMING bound
    EVENT Evacuate KIND other TIMING assumed                     # Chamber evacuated
    RUN Cool DURATION 1800 TEMP 500->25 AMBIENT vacuum PROFILE exponential TAU 300 TIMING bound
    EVENT SE KIND measurement TIMING assumed                     # Spectroscopic ellipsometry at room temperature, two-layer model
}

(100) Ge, 550 °C, 1 atm O2

# (100) Ge, 550 °C, 1 atm O2
SAMPLE Ge100_550 {
    FAMILY: Ge_GeO2
    ORIENTATION: Ge100
    DOPING: n_type
    RESISTIVITY_OHM_CM: 20
    CLEAN: not_stated
}
CONTROL_TRAJECTORY History {
    EVENT Clean KIND clean TIMING assumed                        # Bromine-methanol etch, buffered HF, water rinse
    EVENT Load KIND load TIMING assumed                          # Loaded into the vacuum chamber
    RUN HeatUp DURATION 300 TEMP 25->550 AMBIENT vacuum TIMING bound VARY duration 60..1800
    EVENT Valve KIND measurement TIMING assumed                  # Pump valved off; real-time ellipsometry started
    RUN Fill DURATION 30 TEMP 550->550 AMBIENT mix:O2=0.5 TIMING bound VARY duration 5..300
    RUN Oxidation DURATION 7200 TEMP 550->550 AMBIENT dry_O2 TIMING bound
    EVENT Evacuate KIND other TIMING assumed                     # Chamber evacuated
    RUN Cool DURATION 1800 TEMP 550->25 AMBIENT vacuum PROFILE exponential TAU 300 TIMING bound
    EVENT SE KIND measurement TIMING assumed                     # Spectroscopic ellipsometry at room temperature, two-layer model
}

(100) Ge, 550 °C, 10 min (XPS sample, Table 3)

# (100) Ge, 550 °C, 10 min (XPS sample, Table 3)
SAMPLE Ge100_xps_550C_10min {
    FAMILY: Ge_GeO2
    ORIENTATION: Ge100
    DOPING: n_type
    RESISTIVITY_OHM_CM: 20
    CLEAN: not_stated
}
CONTROL_TRAJECTORY History {
    EVENT Clean KIND clean TIMING assumed                        # Bromine-methanol etch, buffered HF, water rinse
    EVENT Load KIND load TIMING assumed                          # Loaded into the vacuum chamber
    RUN HeatUp DURATION 300 TEMP 25->550 AMBIENT vacuum TIMING bound VARY duration 60..1800
    EVENT Valve KIND measurement TIMING assumed                  # Pump valved off; real-time ellipsometry started
    RUN Fill DURATION 30 TEMP 550->550 AMBIENT mix:O2=0.5 TIMING bound VARY duration 5..300
    RUN Oxidation DURATION 600 TEMP 550->550 AMBIENT dry_O2 TIMING bound
    EVENT Evacuate KIND other TIMING assumed                     # Chamber evacuated
    RUN Cool DURATION 1800 TEMP 550->25 AMBIENT vacuum PROFILE exponential TAU 300 TIMING bound
    EVENT SE KIND measurement TIMING assumed                     # Spectroscopic ellipsometry at room temperature, two-layer model
}

(100) Ge, 400 °C, 540 min (XPS sample, Table 3)

# (100) Ge, 400 °C, 540 min (XPS sample, Table 3)
SAMPLE Ge100_xps_400C_540min {
    FAMILY: Ge_GeO2
    ORIENTATION: Ge100
    DOPING: n_type
    RESISTIVITY_OHM_CM: 20
    CLEAN: not_stated
}
CONTROL_TRAJECTORY History {
    EVENT Clean KIND clean TIMING assumed                        # Bromine-methanol etch, buffered HF, water rinse
    EVENT Load KIND load TIMING assumed                          # Loaded into the vacuum chamber
    RUN HeatUp DURATION 300 TEMP 25->400 AMBIENT vacuum TIMING bound VARY duration 60..1800
    EVENT Valve KIND measurement TIMING assumed                  # Pump valved off; real-time ellipsometry started
    RUN Fill DURATION 30 TEMP 400->400 AMBIENT mix:O2=0.5 TIMING bound VARY duration 5..300
    RUN Oxidation DURATION 32400 TEMP 400->400 AMBIENT dry_O2 TIMING bound
    EVENT Evacuate KIND other TIMING assumed                     # Chamber evacuated
    RUN Cool DURATION 1800 TEMP 400->25 AMBIENT vacuum PROFILE exponential TAU 300 TIMING bound
    EVENT SE KIND measurement TIMING assumed                     # Spectroscopic ellipsometry at room temperature, two-layer model
}

The programs the predictor ran on these histories: scripts/search_growth.py (structure generated from the grammar).

The model itself, scripts/search_growth.py, as it is written. Every equation the simulation uses is in this text; the runtime holds no physics.

#!/usr/bin/env python3
"""
Structure search for oxide growth, across systems. No structure is taken from a paper and nothing is proposed by a language model.
The grammar: the growth rate of the oxide is one of
    linear             a constant
    parabolic          a constant over twice the thickness
    linear_parabolic   B / (2x + B/k), two constants
    power              a constant times thickness to a fitted negative power (0 gives linear, 1 parabolic)
and each constant either has its own activation energy or, for linear_parabolic, the two may share one.
Every candidate is written as an ISL program (kind OxideGrowth), compiled, and run by the runtime.
Two tests, both on measurements the fit does not see:
    A  temperature   SiC, each crystal face: fitted on four temperatures, scored on the held-out 1050 C curve.
    B  time          every growth curve in the library (SiC, planar silicon in steam, a 20% SiGe alloy in steam, a 14%
                     alloy in wet oxygen): constants fitted at that one temperature on all points but the last,
                     scored on the last. This asks which form extrapolates, material by material.
Output: evaluations/growth_search_report.json; the best structure of test A per face is kept under experiments/search_growth/.
"""
import csv
import hashlib
import itertools
import json
import os
import sys
import numpy as np
from scipy.optimize import least_squares
from isl_harness.isl_runtime import ISLRuntime
from scripts import run_song_sic
ROOT = os.path.abspath(os.path.join(os.path.dirname(__file__), ".."))
REPORT = os.path.join(ROOT, "evaluations", "growth_search_report.json")
OUT_DIR = os.path.join(ROOT, "experiments", "search_growth")
T_REF = 1373.15
STRUCTURES = [("linear", "own"), ("parabolic", "own"), ("power", "own"), ("linear_parabolic", "own"), ("linear_parabolic", "shared")]
BOUNDS = {"log_k": (-12.0, -5.0), "E_k": (0.0, 5.0), "log_B": (-17.0, -9.0), "E_B": (0.0, 5.0), "n": (0.0, 2.5), "log_xi": (-8.0, -6.0), "E_xi": (0.0, 3.0)}
START = {"log_k": -8.0, "E_k": 1.5, "log_B": -13.0, "E_B": 1.5, "n": 0.7, "log_xi": -7.0, "E_xi": 0.5}

def names(form, energy, isothermal=False, fit_xi=False): return _names(form, energy, isothermal) + (["log_xi", "E_xi"] if fit_xi == "activated" else ["log_xi"] if fit_xi else [])

def _names(form, energy, isothermal=False): base = {"linear": ["log_k", "E_k"], "parabolic": ["log_B", "E_B"], "power": ["log_k", "E_k", "n"], "linear_parabolic": ["log_B", "E_B", "log_k", "E_k"]}[form] if form == "linear_parabolic" and energy == "shared": base = ["log_B", "E_B", "log_k"] return [n for n in base if not (isothermal and n.startswith("E_"))]

def program(form, energy, sample, ambient, T_C, targets, isothermal=False, dataset="data/derived/song_2004_sic_thickness.csv", sha="0" 64, x_i=2.0e-8, tag="curve", fit_xi=False, start_from_fitted=True, t_ref=T_REF, pressure=False): """fit_xi declares the starting oxide as a fitted parameter; start_from_fitted=False keeps it declared but starts from the given x_i (a measured starting oxide). pressure: the rate is proportional to the oxygen partial pressure (p_O2 / 1 atm), so no oxide grows in nitrogen, argon or vacuum. Every fit was made in oxygen at 1 atm, where the factor is one, so the fitted values are unchanged; proportionality in pressure is the 1965 law's and an assumption for the other forms.""" arr = lambda E: "1" if isothermal else f"exp(-({E} / k_B) (1 / T - 1 / T_ref))" k = f"v_unit exp(ln10 log_k) {arr('E_B' if energy == 'shared' else 'E_k')}" B = f"B_unit exp(ln10 log_B) {arr('E_B')}" rate = {"linear": k, "parabolic": f"{B} / (2 x_ox)", "power": f"{k} exp(-n log(x_ox / x_ref))", "linear_parabolic": f"{B} / (2 x_ox + {B} / ({k}))"}[form] decl = "\n".join(f"PARAMETER {n} {{ VALUE: {START[n]}; UNIT: {'eV' if n.startswith('E_') else '1'}; STATUS: fitted }}" for n in names(form, energy, isothermal, fit_xi)) initial = ("x_unit exp(ln10 log_xi) exp(-(E_xi / k_B) (1 / T_start - 1 / T_ref))" if fit_xi == "activated" else "x_unit exp(ln10 log_xi)" if fit_xi and start_from_fitted else "x_i") tg = "\n".join(f" TARGET: Thickness @ {t:.1f} = {x:.4g} +/- {tol:.4g} abs" for t, x, tol in targets) if pressure: rate = f"(p_O2 / p_ref) ({rate})" return f"""# GENERATED by scripts/search_growth.py. Structure: {form}, activation energy {energy}{', one temperature' if isothermal else ''}. # A candidate from the grammar, not from a paper. Every parameter is fitted. SAMPLE {tag} {{ {sample} }} STATE x_ox {{ UNIT: cm; INITIAL: {initial}; MIN: 0; MAX: 1.0e-3 }} CONTROL T {{ UNIT: K }} {"CONTROL p_O2 { UNIT: atm }" + chr(10) + "PARAMETER p_ref { VALUE: 1.0; UNIT: atm; STATUS: constant }" if pressure else ""} {decl} PARAMETER x_i {{ VALUE: {x_i:.4g}; UNIT: cm; STATUS: constant }} PARAMETER x_unit {{ VALUE: 1.0; UNIT: cm; STATUS: constant }} PARAMETER T_start {{ VALUE: {float(T_C) + 273.15}; UNIT: K; STATUS: constant }} PARAMETER x_ref {{ VALUE: 1.0e-6; UNIT: cm; STATUS: constant }} PARAMETER T_ref {{ VALUE: {t_ref}; UNIT: K; STATUS: constant }} PARAMETER v_unit {{ VALUE: 1.0; UNIT: cm s^-1; STATUS: constant }} PARAMETER B_unit {{ VALUE: 1.0; UNIT: cm^2 s^-1; STATUS: constant }} PARAMETER ln10 {{ VALUE: 2.302585092994046; UNIT: 1; STATUS: constant }} PARAMETER k_B {{ VALUE: 8.617333262e-5; UNIT: eV K^-1; STATUS: constant }} CONTROL_TRAJECTORY History_{tag} {{ RUN Oxidation DURATION {max(t for t, _ in targets):.1f} TEMP {float(T_C):.1f}->{float(T_C):.1f} AMBIENT {ambient} }} OPERATOR F_Growth : OxideGrowth {{ EQUATION: d(x_ox)/dt = {rate} }} OBSERVATION_MODEL Obs {{ TARGET: Thickness; UNIT: cm; EQUATION: Thickness = x_ox }} EVALUATION_CONTRACT Contract_{tag} {{ FROZEN: true DATASET: {dataset} DATASET_DIGEST: sha256:{sha} SPLIT: search {tg} }} """

def fit(form, energy, curves, train_of, isothermal=False, random_starts=6, seed=7, max_nfev=80, fit_xi=False): """curves: list of dicts with rt (parsed program) and targets. train_of(curve) -> indices of targets used to fit.""" ns = names(form, energy, isothermal, fit_xi) lo, hi = [BOUNDS[n][0] for n in ns], [BOUNDS[n][1] for n in ns] def predict(x): ov = dict(zip(ns, [float(v) for v in x])) out = [] for c in curves: try: res = c["rt"].execute_isl_evaluation(parameter_overrides=ov, base_dir=ROOT)["receipt"]["results"] by_time = {round(v["time_s"], 1): v["pred"] for v in res.values()} out.append([by_time[round(t, 1)] for t, _ in c["targets"]]) except Exception: out.append([float("nan")] len(c["targets"])) return out def residuals(x): r = [] for c, p in zip(curves, predict(x)): for i in train_of(c): v = (p[i] - c["targets"][i][1]) / c["targets"][i][2] r.append(v if np.isfinite(v) else 50.0) return r rng, best = np.random.default_rng(seed), None for x0 in [[START[n] for n in ns]] + [[float(rng.uniform(a, b)) for a, b in zip(lo, hi)] for _ in range(random_starts)]: f = least_squares(residuals, x0=x0, bounds=(lo, hi), diff_step=1e-2, max_nfev=max_nfev) if best is None or f.cost < best.cost: best = f return dict(zip(ns, [float(v) for v in best.x])), predict(best.x)

def parse(text): rt = ISLRuntime() rt.parse_isl_source(text) return rt

def test_temperature(job): form, energy, face = job curves_all, sha = run_song_sic.load() curves = [] for (f, T_C), c in sorted(curves_all.items()): if f != face: continue targets = [(t * 3600.0, x, tol) for t, x, tol in sorted(c["points"])] curves.append({"T_C": T_C, "split": c["split"], "targets": targets, "rt": parse(program(form, energy, f"FAMILY: SiC_SiO2; ORIENTATION: 4H_SiC_{face}_face", "dry_O2", T_C, targets, sha=sha, tag=f"SiC_{face}face_{T_C}C"))}) fitted, pred = fit(form, energy, curves, lambda c: range(len(c["targets"])) if c["split"] != "holdout" else []) score = lambda split: [abs(p - t[1]) <= t[2] for c, ps in zip(curves, pred) if (c["split"] == "holdout") == (split == "holdout") for p, t in zip(ps, c["targets"])] return {"structure": f"{form}, energy {energy}", "form": form, "energy": energy, "face": face, "n_fitted": len(fitted), "fitted": fitted, "train_passed": int(sum(score("train"))), "train": len(score("train")), "held_out_passed": int(sum(score("holdout"))), "held_out": len(score("holdout"))}

def test_silicon(job): """Planar silicon, dry: fitted on Massoud at four temperatures; scored on Massoud at 900 C and on another laboratory at 725 and 750 C.""" form, energy = job from scripts import build_silicon_kinetics_dataset as sk rows, sha = sk.build(write=False) curves = [] for data, T_C in sorted({(r["data"], r["T_C"]) for r in rows}): sel = sorted((r for r in rows if r["data"] == data and r["T_C"] == T_C), key=lambda r: r["seconds"]) start = float(sel[0]["x_ox_cm"]) if data == "krzeminski_2007" else None # that paper measures the oxide at time zero targets = [(r["seconds"], float(r["x_ox_cm"]), float(r["tolerance_cm"])) for r in sel if r["seconds"] > 0] curves.append({"data": data, "T_C": T_C, "split": sel[0]["split"], "targets": targets, "rt": parse(program(form, energy, "FAMILY: Si_Planar; ORIENTATION: Si100", "dry_O2", T_C, targets, dataset="data/derived/silicon_dry_kinetics.csv", sha=sha, tag=f"Si_{data}_{T_C}C", fit_xi=True, start_from_fitted=start is None, x_i=start or 1.0e-7))}) fitted, pred = fit(form, energy, curves, lambda c: range(len(c["targets"])) if c["split"] == "train" else [], fit_xi=True) out = {"structure": f"{form}, energy {energy}", "form": form, "energy": energy, "n_fitted": len(fitted), "fitted": fitted} for split in ("train", "holdout", "other_laboratory"): ok = [abs(p - t[1]) <= t[2] for c, ps in zip(curves, pred) if c["split"] == split for p, t in zip(ps, c["targets"])] out[split + "_passed"], out[split] = int(sum(ok)), len(ok) out["train_mean_error_over_tolerance"] = round(float(np.mean([abs(p - t[1]) / t[2] for c, ps in zip(curves, pred) if c["split"] == "train" for p, t in zip(ps, c["targets"])])), 3) out["by_curve"] = [{"data": c["data"], "T_C": c["T_C"], "split": c["split"], "passed": int(sum(abs(p - t[1]) <= t[2] for p, t in zip(ps, c["targets"]))), "of": len(c["targets"]), "mean_signed_error_percent": round(float(np.mean([(p - t[1]) / t[1] for p, t in zip(ps, c["targets"])])) * 100, 1)} for c, ps in zip(curves, pred)] return out

GE_STRUCTURES = [(f, e, start) for f, e in STRUCTURES for start in ("fixed", "activated")]

def test_germanium(job): """Germanium in oxygen at 1 atm: fitted at 400, 450 and 550 C, scored at 500 C. The starting oxide (what the fast first regime leaves) is either one fitted thickness or a thermally activated one.""" form, energy, start = job fit_xi = "activated" if start == "activated" else True from scripts import build_germanium_dataset as gd rows, sha = gd.build(write=False) curves = [] for T_C in sorted({r["T_C"] for r in rows}): sel = sorted((r for r in rows if r["T_C"] == T_C), key=lambda r: r["seconds"]) targets = [(r["seconds"], float(r["x_ox_cm"]), float(r["tolerance_cm"])) for r in sel] curves.append({"T_C": T_C, "split": sel[0]["split"], "targets": targets, "rt": parse(program(form, energy, "FAMILY: Ge_GeO2; ORIENTATION: Ge100", "dry_O2", T_C, targets, dataset="data/derived/germanium_thermal_oxidation.csv", sha=sha, tag=f"Ge_{T_C}C", fit_xi=fit_xi, t_ref=773.15))}) fitted, pred = fit(form, energy, curves, lambda c: range(len(c["targets"])) if c["split"] == "train" else [], fit_xi=fit_xi, random_starts=10) out = {"structure": f"{form}, energy {energy}, starting oxide {start}", "form": form, "energy": energy, "start": start, "n_fitted": len(fitted), "fitted": fitted} for split in ("train", "holdout"): ok = [abs(p - t[1]) <= t[2] for c, ps in zip(curves, pred) if c["split"] == split for p, t in zip(ps, c["targets"])] out[split + "_passed"], out[split] = int(sum(ok)), len(ok) out["train_mean_error_over_tolerance"] = round(float(np.mean([abs(p - t[1]) / t[2] for c, ps in zip(curves, pred) if c["split"] == "train" for p, t in zip(ps, c["targets"])])), 3) out["held_out_points"] = [{"minutes": t[0] / 60, "measured_nm": round(t[1] 1e7, 1), "predicted_nm": round(p 1e7, 1)} for c, ps in zip(curves, pred) if c["split"] == "holdout" for p, t in zip(ps, c["targets"])] return out

def library_curves(): """Every isothermal growth curve with three or more points: (label, sample text, ambient, T_C, [(s, cm, tol)], initial oxide cm).""" out = [] digest = lambda rel: hashlib.sha256(open(os.path.join(ROOT, rel), "rb").read()).hexdigest() song, muller, legoues = "data/derived/song_2004_sic_thickness.csv", "data/derived/muller_2023_series.csv", "data/derived/rabie_2005_legoues_thickness.csv" curves, _ = run_song_sic.load() for (face, T_C), c in sorted(curves.items()): if len(c["points"]) >= 3: out.append((f"SiC {face} face, {T_C} C, dry", f"FAMILY: SiC_SiO2; ORIENTATION: 4H_SiC_{face}_face", "dry_O2", T_C, [(t 3600.0, x, tol) for t, x, tol in sorted(c["points"])], 2.0e-8, song, digest(song))) with open(os.path.join(ROOT, "data", "derived", "muller_2023_series.csv")) as f: rows = [r for r in csv.DictReader(f) if r["shape"] == "planar"] for ge, label, spec in ((0, "planar Si, 850 C, steam", "FAMILY: Si_Planar; SHAPE: planar"), (20, "planar Si0.8Ge0.2, 850 C, steam", "FAMILY: SiGe_oxide; SHAPE: planar; COMPOSITION: Si0.80Ge0.20")): pts = sorted((float(r["seconds"]), float(r["oxide_nm"]) 1e-7, float(r["tolerance_nm"]) 1e-7) for r in rows if int(float(r["germanium_percent"])) == ge) out.append((label, spec, "steam", 850, pts, 1.0e-7, muller, digest(muller))) with open(os.path.join(ROOT, "data", "derived", "rabie_2005_legoues_thickness.csv")) as f: rows = list(csv.DictReader(f)) tk = next(k for k in rows[0] if "time" in k.lower()); xk = next(k for k in rows[0] if "cm" in k.lower() and "tol" not in k.lower()); ok = next((k for k in rows[0] if "tol" in k.lower()), None) scale = 60.0 if "min" in tk.lower() else 3600.0 if "h" in tk.lower() else 1.0 out.append(("planar Si0.86Ge0.14, 800 C, wet oxygen", "FAMILY: SiGe_oxide; SHAPE: planar; COMPOSITION: Si0.86Ge0.14", "wet_O2", 800, sorted((float(r[tk]) scale, float(r[xk]), float(r[ok]) if ok else 0.1 * float(r[xk])) for r in rows), 1.0e-7, legoues, digest(legoues))) return out

def test_time(job): form, index = job label, spec, ambient, T_C, targets, x_i, dataset, sha = library_curves()[index] curve = {"targets": targets, "rt": parse(program(form, "own", spec, ambient, T_C, targets, isothermal=True, x_i=x_i, tag=f"curve_{index}", dataset=dataset, sha=sha))} fitted, pred = fit(form, "own", [curve], lambda c: range(len(c["targets"]) - 1), isothermal=True, random_starts=4, max_nfev=60) last = targets[-1] return {"curve": label, "form": form, "points_fitted": len(targets) - 1, "n_fitted": len(fitted), "fitted": fitted, "last_measured_nm": last[1] 1e7, "last_predicted_nm": float(pred[0][-1]) 1e7, "passed": bool(abs(pred[0][-1] - last[1]) <= last[2]), "determined": len(targets) - 1 > len(fitted)}

def where_the_forms_part(by_face, sha): """Thickness the power law and the linear-parabolic law predict, each as fitted, at times beyond the SiC data (which end near 6 h).""" rows = [] for face in ("C", "Si", "a"): for T_C in (1050, 1150): row = {"face": face, "T_C": T_C} for hours in (6, 24, 96): pred = {} for r in by_face: if r["face"] == face and r["structure"] in ("power, energy own", "linear_parabolic, energy own"): rt = parse(program(r["form"], r["energy"], f"FAMILY: SiC_SiO2; ORIENTATION: 4H_SiC_{face}_face", "dry_O2", T_C, [(hours 3600.0, 1e-5, 1e-6)], sha=sha, tag="design")) pred[r["form"]] = next(iter(rt.execute_isl_evaluation(parameter_overrides=r["fitted"], base_dir=ROOT)["receipt"]["results"].values()))["pred"] 1e7 row[f"{hours}h"] = {"power_nm": round(pred["power"], 1), "linear_parabolic_nm": round(pred["linear_parabolic"], 1), "differ_percent": round(100 * abs(pred["power"] - pred["linear_parabolic"]) / pred["linear_parabolic"], 1)} rows.append(row) return {"what": "Thickness predicted by the power law and by the linear-parabolic law, each as fitted, at times beyond the data (the data end near 6 h). A contract allows 10%.", "rows": rows}

def run(workers=1): import multiprocessing jobs_a = [(form, energy, face) for form, energy in STRUCTURES for face in ("C", "Si", "a")] jobs_b = [(form, i) for form in ("linear", "parabolic", "power", "linear_parabolic") for i in range(len(library_curves()))] if workers > 1: with multiprocessing.Pool(workers) as pool: a, b, si, ge = pool.map(test_temperature, jobs_a), pool.map(test_time, jobs_b), pool.map(test_silicon, STRUCTURES), pool.map(test_germanium, GE_STRUCTURES) else: a, b, si, ge = [test_temperature(j) for j in jobs_a], [test_time(j) for j in jobs_b], [test_silicon(j) for j in STRUCTURES], [test_germanium(j) for j in GE_STRUCTURES] # Selection uses the fitted temperatures only (most training targets passed, then fewest parameters). Held-out temperatures and other # laboratories are scored afterwards and are not used to choose, so they stay an assessment of the structure that was chosen. key = lambda r: (-r["train_passed"], r["n_fitted"], r["train_mean_error_over_tolerance"]) select = lambda rows: min(rows, key=key) RULE = "most training targets passed, then fewest parameters, then smallest training error; held-out scores are not used to choose" ge.sort(key=key) si.sort(key=key) by_structure = {} for r in a: s = by_structure.setdefault(r["structure"], {"structure": r["structure"], "n_fitted_per_face": r["n_fitted"], "held_out_passed": 0, "held_out": 0, "train_passed": 0, "train": 0}) for k in ("held_out_passed", "held_out", "train_passed", "train"): s[k] += r[k] ranking_a = sorted(by_structure.values(), key=lambda s: (-s["train_passed"], s["n_fitted_per_face"])) by_form = {} for r in b: f = by_form.setdefault(r["form"], {"form": r["form"], "passed": 0, "of": 0, "passed_where_determined": 0, "determined": 0, "by_material": {}}) f["passed"] += r["passed"]; f["of"] += 1 f["determined"] += r["determined"]; f["passed_where_determined"] += r["passed"] and r["determined"] m = "SiC" if r["curve"].startswith("SiC") else r["curve"] f["by_material"].setdefault(m, [0, 0]); f["by_material"][m][0] += r["passed"]; f["by_material"][m][1] += 1 ranking_b = sorted(by_form.values(), key=lambda f: -f["passed"]) os.makedirs(OUT_DIR, exist_ok=True) for old in os.listdir(OUT_DIR): os.remove(os.path.join(OUT_DIR, old)) best = ranking_a[0]["structure"] curves_all, sha = run_song_sic.load() for r in a: if r["structure"] == best and r["held_out"]: T_C = next(T for (f, T), c in curves_all.items() if f == r["face"] and c["split"] == "holdout") c = curves_all[(r["face"], T_C)] text = program(r["form"], r["energy"], f"FAMILY: SiC_SiO2; ORIENTATION: 4H_SiC_{r['face']}_face", "dry_O2", T_C, [(t * 3600.0, x, tol) for t, x, tol in sorted(c["points"])], sha=sha, tag=f"SiC_{r['face']}face_{T_C}C") for n, v in r["fitted"].items(): text = text.replace(f"PARAMETER {n} {{ VALUE: {START[n]};", f"PARAMETER {n} {{ VALUE: {v:.6g};") with open(os.path.join(OUT_DIR, f"best_{r['form']}_{r['energy']}_SiC_{r['face']}face.isl"), "w") as f: f.write(text) report = {"what": "Structure search for oxide growth over ISL programs generated from a grammar. No structure is taken from a paper.", "data_verification": "Agent-read values. Unverified, except Mueller's, which were extracted from the vector drawing.", "grammar": {"forms": ["linear", "parabolic", "power", "linear_parabolic"], "activation_energy": ["own", "shared"]}, "test_A_temperature_SiC": {"what": "fitted on the other temperatures of each face, scored on the held-out 1050 C curves", "ranking": ranking_a, "by_face": a}, "where_the_two_surviving_forms_part": where_the_forms_part(a, sha), "test_C_planar_silicon_dry": {"what": "fitted on Massoud's data at 800, 850, 950 and 1000 C (as redrawn by Krzeminski et al.); scored on Massoud at 900 C and on " "Krzeminski's own measurements at 725 and 750 C. The starting oxide is one further fitted number.", "selection_rule": RULE, "selected": select(si), "ranking": si}, "test_D_germanium": {"what": "germanium in oxygen at 1 atm (Wang, Hu and Irene 1994, read from a scanned figure): fitted at 400, 450 and 550 C, scored at 500 C. " "Rate constants are referred to 500 C.", "selection_rule": RULE, "selected": select(ge), "ranking": ge}, "test_B_time_all_materials": {"what": "fitted at one temperature on all points but the last of each curve, scored on the last", "ranking": ranking_b, "by_curve": b}} with open(REPORT, "w") as f: json.dump(report, f, indent=2) return report

if __name__ == "__main__": workers = int(sys.argv[sys.argv.index("--workers") + 1]) if "--workers" in sys.argv else 1 rep = run(workers) print("A. SiC, held-out temperature:") for s in rep["test_A_temperature_SiC"]["ranking"]: print(f" {s['structure']:34s} held-out {s['held_out_passed']}/{s['held_out']} train {s['train_passed']}/{s['train']} parameters per face {s['n_fitted_per_face']}") print("C. planar silicon, dry:") for s in rep["test_C_planar_silicon_dry"]["ranking"]: print(f" {s['structure']:34s} train {s['train_passed']}/{s['train']} 900 C held out {s['holdout_passed']}/{s['holdout']} other laboratory {s['other_laboratory_passed']}/{s['other_laboratory']} parameters {s['n_fitted']}") print("D. germanium:") for s in rep["test_D_germanium"]["ranking"]: print(f" {s['structure']:58s} train {s['train_passed']}/{s['train']} 500 C held out {s['holdout_passed']}/{s['holdout']} parameters {s['n_fitted']} {s['held_out_points']}") print("B. last point of every curve:") for f in rep["test_B_time_all_materials"]["ranking"]: print(f" {f['form']:18s} {f['passed']}/{f['of']} {f['by_material']}") ```

Our simulation against the paper

8 within, 2 near, 3 off, 3 not modelled of 16 numbers. Tolerance is stated per row; nothing was fitted for this page.

Aligning the simulation: what the unstated things would have to be

Every step whose timing or temperature the paper leaves open carries a range in the history above (VARY). scripts/align_study.py searched those ranges against the paper's numbers through the same predictor, first one unknown at a time, then all together. No parameter a program declares fitted was touched. The unknowns were chosen on 400C, 450C, 550C only; every other number is a test they played no part in.

Together: HeatUp duration 300 s (was 300); Fill duration 30 s (was 30).

As reproduced: 8 of 16 within tolerance, median error 2.7%. Aligned: 8 of 16 within tolerance, median error 2.7%. On the held-out samples only: 0 of 6 within tolerance, median error 22.0% as reproduced, 0 of 6 within tolerance, median error 22.0% aligned.

Choices of program or material property, each scored on its own and on top of the aligned history:

Two unknowns were searched: the heat-up in vacuum (1 to 30 min) and the fill of the closed chamber with oxygen (5 s to 5 min), chosen on the 400, 450 and 550 °C curves and tested on 500 °C and the two Table 3 samples. Neither moves anything that matters. All 8 points of the three fitted curves stay within tolerance over both whole ranges (median error 1.3 to 1.5% for the heat-up, 1.3 to 2.4% for the fill); the search picks the shortest of each because the program grows a few hundredths of a nanometre in vacuum and at the fill's half pressure as if in full oxygen. The held-out numbers are unchanged: 0 of 5 within tolerance, median error 22%, before and after.

So the data cannot tell us the heat-up time or the fill time, and the program could not use them if they could: it reads neither the gas nor the pressure. The variant with the heat-up in oxygen, standing for a conventional furnace, gives the same thicknesses to three figures, but every number is then labelled outside the tested scope, because the predictor takes the start temperature of the first oxygen step (25 °C) as the oxidation temperature.

The starting oxide, the third unknown the paper leaves open (about 0 nm in Figure 1a, about 1 nm in Figure 1b), could not be searched: it is a fitted parameter inside the program, not a property of the sample.

2 near, 3 off, 1 not modelled of 6 numbers with the aligned settings (the samples used to choose the unknowns are left out).

What alignment cannot fix. At the held-out 500 °C the program is 11, 15 and 22% low at 38, 124 and 597 min, and the error grows with time: the power law fitted on 400, 450 and 550 °C does not carry its temperature dependence to 500 °C. That is a model question, not a timing one. The two Table 3 samples are 74% and 49% above the paper's printed thicknesses (8 nm at 550 °C, 10 min; 3 nm at 400 °C, 540 min), but the paper's own Figure 1 is also above them (about 10 nm and 4.4 nm). They were probably measured with another optical model or are other samples; the table does not say.

The aligned score is not a test on the fitted temperatures: the program's four numbers were fitted there. Only the 500 °C curve and the Table 3 rows are held out.

What our language and ontology cannot yet say

• A starting oxide on the sample. The paper's two figures disagree on it (0 or 1 nm) and our program holds it as a fitted parameter. A SAMPLE … STARTING_OXIDE <nm> field that growth programs read would make it an unknown the alignment can search and a fact a paper can state.
• Growth programs that read the oxygen pressure. The germanium program grows oxide at full rate in vacuum and at half an atmosphere. Every growth law should multiply by a pressure dependence (linear in p_O2 for the parabolic term at least); then the heat-up, the fill and the cooling in vacuum become testable instead of inert.
• A pressure that changes within a step. “Oxygen was admitted” into a closed chamber at 150 sccm is a ramp of pressure, as a ramp of temperature is. RUN … PRESSURE 0->1 (or a pressure PROFILE) would let the fill be written as stated, with the chamber volume as the unknown.
• A measured clock origin. The paper's time zero is the admission of oxygen, which is the start of the fill, not of our Oxidation step. A curve should name the event its clock starts at (CLOCK FROM Fill), so a fill or ramp in front does not shift every point.
• Two observables for one oxide. Total thickness of a two-layer optical model, the GeO2 layer alone, the interface layer (GeO2 + amorphous Ge, 1.4 to 4.8 nm), and a single-layer real-time value are four numbers for one sample here. The ontology has one Thickness. An interface-layer observable and an observation model per technique would let Table 1 be scored.
• Scope from the oxidation, not from the first oxygen step. The predictor labels a history by the start temperature of its first oxygen step, so a ramp in oxygen makes a 550 °C oxidation “outside the tested scope”.
• A germanium clean. The vocabulary of pre-oxidation cleans is silicon's; bromine-methanol with buffered HF cannot be named.
• The program as a file. The selected germanium law is generated in memory by scripts/search_growth.py; writing it to experiments/ as the SiC ones are would let this page print the model it ran.

With heating and cooling transitions

No paper we hold states how its wafers heated up or cooled down. The project owner specified scenarios for one Vitkavage wafer (isl_harness/vitkavage_transitions.py); isl_harness/thermal_transitions.py applies the same numbers to every sample here: a heat-up in nitrogen before the first furnace step (a hot push of 1 to 5 min, or a cold-furnace ramp of 60 to 120 min), the owner's radiative cooling knots after the last furnace step or after the paper's own pull (fast, or every time multiplied by 5), and ramps around an anneal at 300 to 500 °C. The paper's own steps are kept as written. All transitions are assumptions, written in the history with TIMING assumed.

The paper's own conclusions

What we would like Emily to check

1. How long would a chamber of this kind take to reach 1 atm at 150 sccm, and does the thickness at a few minutes then mean anything?
2. Does “Figure 1d … conventional oxidation furnace” refer to the long-time thermal panel (1b)?
3. Which thickness should a model be scored on: the total of the two layers, or the GeO2 layer alone?
4. Table 3 and Figure 1 disagree for the same conditions. Are the XPS samples different samples, measured with a single-layer model?