Rewrite in Rust with Python/Jupyter frontend, rename to NEGForge

[reneotten]

Summary

• Ports the MATLAB 1D ballistic-MOSFET / NEGF simulator (from the testing branch, the furthest-ahead branch with the actual physics code — master only had early drafts) to a Rust engine exposed to Python via PyO3, with a Jupyter notebook frontend.
• crates/negforge-core: electrostatics (tridiagonal Poisson-like solve), ballistic Landauer current, and an O(N) recursive-Green's-function NEGF implementation (replacing the original's O(N³) dense inversion), validated against dense reference solvers in tests.
• New: a self-consistent Poisson↔NEGF feedback loop the original MATLAB code never implemented, plus two bug fixes in the NEGF charge formula needed to make that loop physically meaningful (energy-dependent Fermi weighting; per-contact band-edge masking). See the README's "Deviations from the original" section for full rationale, including the charge-density unit-scaling fix required to keep the iteration numerically stable.
• crates/negforge-py: PyO3 bindings; python/negforge: a thin numpy-friendly wrapper package.
• notebooks/negforge_demo.ipynb: interactive frontend reproducing all of the original plots (potential profile, I-Vg/I-Vds curves, subthreshold swing, NEGF local-density-of-states heatmap) plus the new self-consistent solve, executed end-to-end to confirm it works.
• legacy_matlab/: original MATLAB code kept for reference.
• Project renamed from quantumsim to NEGForge (crates, Python package, notebook, docs) per follow-up request.

Test plan

• cargo test --workspace — 16 tests pass (unit tests per module, cross-checked against dense reference solvers; an integration test at realistic device scale confirming the self-consistent loop converges and current increases with gate bias)
• cargo clippy --workspace --all-targets — clean
• maturin develop --release builds and installs the negforge Python package
• notebooks/negforge_demo.ipynb executed end-to-end via jupyter nbconvert --execute with no errors; physics sanity-checked (subthreshold swing ≈ 66 mV/decade, close to the 60 mV/decade ideal thermal limit at 300 K)

---

Generated by Claude Code

[Copilot]

Pull request overview

Ports the legacy MATLAB 1D ballistic-MOSFET / NEGF simulator to a Rust core engine with PyO3 bindings and a Jupyter-based Python frontend, while renaming the project to NEGForge and adding a self-consistent Poisson↔NEGF feedback loop.

Changes:

• Introduces negforge-core (Rust) implementing electrostatics, NEGF, charge density, and a self-consistent Poisson↔NEGF solve with tests.
• Adds negforge-py + python/negforge wrapper package, plus a demo notebook for interactive plotting/sweeps.
• Preserves the original MATLAB sources under legacy_matlab/ for reference and documents deliberate deviations/bug fixes in the new README.

Reviewed changes

Copilot reviewed 24 out of 27 changed files in this pull request and generated 5 comments.

Show a summary per file

File	Description
README.md	New top-level documentation: architecture, build/run steps, physics model, and deviations/bug fixes.
python/negforge/init.py	Python wrapper API (Device, IVCurve) around the PyO3 extension module.
pyproject.toml	Python packaging configuration using maturin and numpy dependency.
notebooks/negforge_demo.ipynb	End-to-end demo notebook reproducing plots and showcasing self-consistent solve.
legacy_matlab/README.md	Explains purpose/limitations of preserved MATLAB reference code.
legacy_matlab/quantumsim.m	Preserved “most complete” legacy MATLAB class implementation for reference.
legacy_matlab/simulation.m	Preserved procedural MATLAB variant (reference only).
legacy_matlab/simulation_v2.m	Preserved procedural MATLAB sweep variant (reference only).
legacy_matlab/simulation_v3.m	Preserved procedural MATLAB sweep variant (reference only).
legacy_matlab/untitled3.m	Preserved legacy MATLAB script snippet (reference only).
legacy_matlab/aufgabe1e.m	Preserved simpler Poisson-only exercise (reference only).
crates/negforge-core/Cargo.toml	Defines the Rust core crate and its dependencies/test deps.
crates/negforge-core/src/lib.rs	Exposes core modules and public re-exports for Rust API consumers.
crates/negforge-core/src/constants.rs	Adds physical constants module replacing missing util.const.
crates/negforge-core/src/error.rs	Defines core error types (notably non-convergence / invalid params).
crates/negforge-core/src/tridiag.rs	Implements tridiagonal solvers + reference tests.
crates/negforge-core/src/device.rs	Implements device params, electrostatics, and ballistic current model.
crates/negforge-core/src/negf.rs	Implements O(N) NEGF sweep (diag + boundary columns) + reference tests.
crates/negforge-core/src/charge.rs	Implements NEGF-derived electron density with documented bug fixes.
crates/negforge-core/src/selfconsistent.rs	Adds new self-consistent Poisson↔NEGF loop and tests.
crates/negforge-core/src/sweep.rs	Adds bias sweeps and subthreshold-swing extraction utilities + tests.
crates/negforge-core/tests/self_consistent_realistic_device.rs	Integration test verifying convergence and current trend at default device scale.
crates/negforge-py/Cargo.toml	Defines the PyO3 bindings crate as a cdylib.
crates/negforge-py/src/lib.rs	PyO3 bindings exposing the Rust Device to Python and sweep/LDOS helpers.
Cargo.toml	Introduces workspace layout and release profile settings.
Cargo.lock	Locks Rust dependencies for reproducible builds.
.gitignore	Ignores Rust build artifacts, Python env/build outputs, and notebook checkpoints.

Comments suppressed due to low confidence (6)

python/negforge/init.py:97

• After set_v_g updates the underlying Rust device, the electrostatic potential is not recomputed. This leaves psi_f/psi_0 stale, so calc_current() and local_density_of_states() can silently use the previous bias point unless the caller remembers to call calc_potential() manually.

    def set_v_g(self, v: float) -> "Device":
        self._inner.set_v_g(v)
        return self

python/negforge/init.py:101

• After set_l_ch changes device geometry, the electrostatic potential is not recomputed. This leaves psi_f/psi_0 inconsistent with the new grid/geometry, so subsequent calc_current() calls can be wrong unless calc_potential() is invoked manually.

    def set_l_ch(self, l: float) -> "Device":
        self._inner.set_l_ch(l)
        return self

crates/negforge-core/src/sweep.rs:102

• subthreshold_swing takes log10() of the raw current values. If any current in the fit window is non-positive (or non-finite), log10() yields NaN/-inf and the fit result becomes NaN (often without an obvious error). Filtering to finite, positive current points (and returning NaN when insufficient data remains) makes this utility more robust.

pub fn subthreshold_swing(points: &[IvPoint], fit_v_min: f64, fit_v_max: f64) -> f64 {
    let (xs, ys): (Vec<f64>, Vec<f64>) = points
        .iter()
        .filter(|p| p.voltage >= fit_v_min && p.voltage <= fit_v_max)
        .map(|p| (p.voltage, p.current.log10()))
        .unzip();
    let (slope, _intercept) = linear_fit(&xs, &ys);
    1.0 / slope
}

crates/negforge-core/src/device.rs:261

• set_v_g resets the device drive terms via init_vectors() but does not recompute psi_f/psi_0. Since calc_current() uses psi_0 as the lower integration bound, calling calc_current() immediately after set_v_g() can silently use a stale potential from the previous bias point.

    pub fn set_v_g(&mut self, v: f64) {
        self.params.v_g = v;
        self.init_vectors();
    }

crates/negforge-core/src/device.rs:270

• set_l_ch rebuilds the grid and resets vectors, but does not recompute psi_f/psi_0. This can leave the device in an inconsistent state where subsequent calc_current() uses an integration lower bound computed on the old geometry.

    pub fn set_l_ch(&mut self, l: f64) {
        self.params.l_ch = l;
        let (n, n_left, n_right) = Self::compute_grid(&self.params);
        self.n = n;
        self.n_left = n_left;
        self.n_right = n_right;
        self.init_vectors();
    }

crates/negforge-core/src/tridiag.rs:67

• solve_complex panics for n == 0 (underflow in n - 1) and for n == 1 (c_prime[0] out of bounds). Handle the 1x1 case explicitly (and assert n > 0) to avoid unexpected panics when solving degenerate systems.

    let n = diag.len();
    assert_eq!(sub.len(), n - 1);
    assert_eq!(sup.len(), n - 1);
    assert_eq!(rhs.len(), n);

    let mut c_prime = vec![Complex64::new(0.0, 0.0); n - 1];
    let mut d_prime = vec![Complex64::new(0.0, 0.0); n];

    c_prime[0] = Complex64::new(sup[0], 0.0) / diag[0];
    d_prime[0] = rhs[0] / diag[0];

---

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.