Adom Fields

@@ -0,0 +1,90 @@+//! Plated-hole geometry and distributed layer coupling. No conducting disc in the drill.+use crate::board::Via;+use std::collections::BTreeMap;+const PLATING_MM: f64 = 0.025;++pub fn transform(v: &Via, x: f64, y: f64) -> (f64, f64) {+    let (s, c) = v.drill_angle.to_radians().sin_cos();+    (v.x + x * c + y * s, v.y - x * s + y * c)+}+pub fn axis(v: &Via) -> ((f64, f64), (f64, f64)) {+    let a = (v.drill_long.unwrap_or(v.drill) - v.drill) / 2.0;+    (transform(v, -a, 0.0), transform(v, a, 0.0))+}+/// Uniform arc-length quadrature at the middle of the plated wall. The weights+/// sum to one; multiplying by sigma * barrel area / length preserves conductance.+pub fn weights(v: &Via, x0: f64, y0: f64, pitch: f64, w: usize, h: usize) -> Vec<(usize, f64)> {+    let a = (v.drill_long.unwrap_or(v.drill) - v.drill) / 2.0;+    let r = v.drill / 2.0 + PLATING_MM / 2.0;+    let arc = std::f64::consts::PI * r;+    let perimeter = 4.0 * a + 2.0 * arc;+    let samples = ((perimeter / pitch * 32.0).ceil() as usize).max(128);+    let mut bins = BTreeMap::<usize, usize>::new();+    for j in 0..samples {+        let t = (j as f64 + 0.5) * perimeter / samples as f64;+        let (x, y) = if t < 2.0 * a {+            (-a + t, -r)+        } else if t < 2.0 * a + arc {+            let q = (t - 2.0 * a) / r - std::f64::consts::FRAC_PI_2;+            (a + r * q.cos(), r * q.sin())+        } else if t < 4.0 * a + arc {+            (a - (t - 2.0 * a - arc), r)+        } else {+            let q = (t - 4.0 * a - arc) / r + std::f64::consts::FRAC_PI_2;+            (-a + r * q.cos(), r * q.sin())+        };+        let (x, y) = transform(v, x, y);+        let c = ((x - x0) / pitch).floor() as isize;+        let row = ((y - y0) / pitch).floor() as isize;+        if c >= 0 && row >= 0 && c < (w as isize) && row < (h as isize) {+            *bins.entry(row as usize * w + c as usize).or_default() += 1;+        }+    }+    bins.into_iter()+        .map(|(i, n)| (i, n as f64 / samples as f64))+        .collect()+}+#[cfg(test)]+mod tests {+    use super::*;+    fn via() -> Via {+        Via {+            net: "N".into(),+            x: 3.0,+            y: 3.0,+            size: 1.2,+            drill: 0.6,+            drill_long: None,+            drill_angle: 0.0,+        }+    }+    #[test]+    fn coupling_is_conserved_and_avoids_center_at_both_meshes() {+        for pitch in [0.2, 0.1, 0.05] {+            let n = (6.0 / pitch) as usize;+            let ws = weights(&via(), 0.0, 0.0, pitch, n, n);+            assert!((ws.iter().map(|x| x.1).sum::<f64>() - 1.0).abs() < 1e-12);+            let center = (3.0 / pitch) as usize * n + (3.0 / pitch) as usize;+            assert!(!ws.iter().any(|(i, _)| *i == center));+            assert!(ws.len() > 4);+        }+    }+    #[test]+    fn rotated_slot_preserves_conductance_and_changes_axis() {+        let mut v = via();+        v.drill_long = Some(1.6);+        v.drill_angle = 90.0;+        let (a, b) = axis(&v);+        assert!((a.0 - b.0).abs() < 1e-12);+        assert!((a.1 - b.1).abs() > 0.99);+        assert!(+            (weights(&v, 0.0, 0.0, 0.1, 60, 60)+                .iter()+                .map(|x| x.1)+                .sum::<f64>()+                - 1.0)+                .abs()+                < 1e-12+        );+    }+}

@@ -9,7 +9,11 @@ pub struct Poly { pub net: String, pub layer: String, pub pts: Vec<(f64, f64)> } #[derive(Debug, Clone)] pub struct Track { pub net: String, pub layer: String, pub a: (f64, f64), pub b: (f64, f64), pub w: f64 } #[derive(Debug, Clone)]-pub struct Via { pub net: String, pub x: f64, pub y: f64, pub size: f64, pub drill: f64 }+pub struct Via { pub net: String, pub x: f64, pub y: f64, pub size: f64, pub drill: f64, /** a slotted hole's long dimension */ pub drill_long: Option<f64>, /** the slot's angle, degrees */ pub drill_angle: f64 }+impl Via {+    /// The plated barrel's copper cross-section: a circular or slotted hole, 25 µm plating (PR 5).+    pub fn barrel_area_mm2(&self) -> f64 { let t = 0.025; let per = std::f64::consts::PI * self.drill + 2.0 * (self.drill_long.unwrap_or(self.drill) - self.drill); per * t + std::f64::consts::PI * t * t }+} #[derive(Debug, Clone)] pub struct Pad { pub reference: String, pub pad: String, pub net: String, pub x: f64, pub y: f64, pub layers: Vec<String>, pub area: f64, /** the symbol pin's name when the file has it (KiCad pinfunction: VIN, VOUT_3, GND) */ pub func: String } @@ -30,6 +34,8 @@ pub struct Board {     pub pad_polys: Vec<Poly>,     pub arc_nets: Vec<String>,     pub layer_names: Vec<String>,+    /// plated through-hole pads' barrels, joining the copper layers like vias (PR 5)+    pub plated_holes: Vec<Via>,     /// pours to fill on the solve grid (Fusion/EAGLE boards; empty for KiCad, whose fill is in the file)     pub pours: Vec<Pour>,     /// pours whose fill the file DOES carry as outline-with-holes islands (Fusion's polygonfilldetails)@@ -108,8 +114,11 @@ impl Board {             b.tracks.push(Track { net: s.kid("net").map(&net_name).unwrap_or_default(), layer: s.kid("layer").and_then(|l| l.atom(1)).unwrap_or("").to_string(), a: (st.num(1).unwrap_or(0.0), st.num(2).unwrap_or(0.0)), b: (en.num(1).unwrap_or(0.0), en.num(2).unwrap_or(0.0)), w: s.kid("width").and_then(|w| w.num(1)).unwrap_or(0.2) });         }         for v in root.kids("via") {+            // the field model joins every copper layer through a via: a blind, buried or micro via is refused (PR 5)+            if v.atom(1).is_some_and(|k| k == "blind" || k == "micro") { return Err("blind, buried or micro vias are not supported: the field model joins every copper layer through a via".into()); }+            if let Some(span) = v.kid("layers") { if span.atom(1).is_some_and(|a| a != "F.Cu") || span.atom(2).is_some_and(|z| z != "B.Cu") { return Err("only through vias (F.Cu to B.Cu) are supported".into()); } }             let Some(at) = v.kid("at") else { continue };-            b.vias.push(Via { net: v.kid("net").map(&net_name).unwrap_or_default(), x: at.num(1).unwrap_or(0.0), y: at.num(2).unwrap_or(0.0), size: v.kid("size").and_then(|s| s.num(1)).unwrap_or(0.6), drill: v.kid("drill").and_then(|d| d.num(1)).unwrap_or(0.3) });+            b.vias.push(Via { net: v.kid("net").map(&net_name).unwrap_or_default(), x: at.num(1).unwrap_or(0.0), y: at.num(2).unwrap_or(0.0), size: v.kid("size").and_then(|s| s.num(1)).unwrap_or(0.6), drill: v.kid("drill").and_then(|d| d.num(1)).unwrap_or(0.3), drill_long: None, drill_angle: 0.0 });         }         for fp in root.kids("footprint") {             let reference = fp.kids("property").iter().find(|p| p.atom(1) == Some("Reference")).and_then(|p| p.atom(2)).unwrap_or("?").to_string();@@ -129,6 +138,16 @@ impl Board {                 let (sw, sh) = (p.kid("size").and_then(|s| s.num(1)).unwrap_or(1.0), p.kid("size").and_then(|s| s.num(2)).unwrap_or(1.0));                 let kind = p.atom(2).unwrap_or("smd");                 let shape = p.atom(3).unwrap_or("rect");+                if kind == "thru_hole" {+                    if let Some(dn) = p.kid("drill") {+                        let (dr, dl) = if dn.atom(1) == Some("oval") { let a = dn.num(2).unwrap_or(0.0); let z = dn.num(3).unwrap_or(a); (a.min(z), Some(a.max(z))) } else { (dn.num(1).unwrap_or(0.0), None) };+                        // a slot along its own axis (PR 9); an offset drill moves the hole with the pad's rotation+                        let ang = pa + if dn.atom(1) == Some("oval") && dn.num(3) > dn.num(2) { 90.0 } else { 0.0 };+                        let (ox, oy) = dn.kid("offset").map(|o| (o.num(1).unwrap_or(0.0), o.num(2).unwrap_or(0.0))).unwrap_or((0.0, 0.0));+                        let (qx, qy) = rot(ox, -oy, pa);+                        if dr > 0.0 { b.plated_holes.push(Via { net: net.clone(), x: x + qx, y: y - qy, size: sw.min(sh), drill: dr, drill_long: dl, drill_angle: ang }); }+                    }+                }                 let lay: Vec<String> = p.kid("layers").map(|l| l.list().iter().skip(1).filter_map(|x| if let Sx::Atom(a) = x { Some(a.clone()) } else { None }).collect()).unwrap_or_default();                 // only onto copper the board has: a plated pad on a two-layer board is on F.Cu and B.Cu, never In1/In2                 let real = |l: &&&str| b.copper.is_empty() || b.copper.iter().any(|c| c == **l);@@ -226,3 +245,33 @@ pub fn check_copper(copper: &[String], declared: &[String]) -> Result<(), String         assert_eq!(s.len(), 5); assert!((s[0] - 0.2).abs() < 1e-12 && (s[1] - 1.065).abs() < 1e-12 && (s[2] - 0.07).abs() < 1e-12 && (s[4] - 0.07).abs() < 1e-12);     } }++#[cfg(test)] mod barrel_tests {+    use super::*;+    // PR 5: a two-layer board's plated pad and via stay on F.Cu and B.Cu, and the pad's barrel is kept+    #[test] fn two_layer_plated_pad_and_via() {+        let p = std::env::temp_dir().join(format!("fields-barrel-{}.kicad_pcb", std::process::id()));+        std::fs::write(&p, r#"(kicad_pcb (layers (0 "F.Cu" signal) (31 "B.Cu" signal) (44 "Edge.Cuts" user))+          (net 1 "GND") (gr_rect (start 0 0) (end 10 10) (layer "Edge.Cuts"))+          (via (at 5 5) (size 1) (drill 0.3) (layers "F.Cu" "B.Cu") (net 1))+          (footprint "test" (property "Reference" "J1") (at 3 3)+            (pad "1" thru_hole circle (at 0 0) (size 2 2) (drill 1) (layers "*.Cu" "*.Mask") (net 1 "GND"))))"#).unwrap();+        let b = Board::load(p.to_str().unwrap()).unwrap(); let _ = std::fs::remove_file(&p);+        assert_eq!(b.pads[0].layers, vec!["F.Cu".to_string(), "B.Cu".to_string()]);+        assert_eq!(b.plated_holes.len(), 1); assert_eq!(b.plated_holes[0].drill, 1.0);+        let mut ras = crate::raster::Raster::new(&b, 0.1); ras.draw_board(&b);+        assert!(ras.masks.keys().all(|(_, l)| l == "F.Cu" || l == "B.Cu"));+    }+    #[test] fn circular_and_slotted_barrels() {+        let round = Via { net: "GND".into(), x: 0., y: 0., size: 2., drill: 1., drill_long: None, drill_angle: 0.0 };+        assert!((round.barrel_area_mm2() - 0.08050331175).abs() < 1e-9);+        let slot = Via { drill: 0.6, drill_long: Some(1.7), ..round };+        assert!((slot.barrel_area_mm2() - 0.10408738521).abs() < 1e-9);+    }+    #[test] fn blind_via_refused() {+        let p = std::env::temp_dir().join(format!("fields-blind-{}.kicad_pcb", std::process::id()));+        std::fs::write(&p, r#"(kicad_pcb (layers (0 "F.Cu" signal) (1 "In1.Cu" signal) (2 "In2.Cu" signal) (31 "B.Cu" signal))+          (gr_rect (start 0 0) (end 10 10) (layer "Edge.Cuts")) (via blind (at 5 5) (size 0.6) (drill 0.3) (layers "F.Cu" "In1.Cu") (net 1)))"#).unwrap();+        let r = Board::load(p.to_str().unwrap()); let _ = std::fs::remove_file(&p); assert!(r.is_err());+    }+}

@@ -24,6 +24,8 @@ pub struct Spec {     #[serde(default)] pub copperUmByLayer: BTreeMap<String,f64>,     #[serde(default)] pub signals: crate::signals::Config,     #[serde(default)] pub heatsinks: Vec<crate::heatsink::Block>,+    /// a power device cooled off the board: junction, case and sink nodes with supplied resistances (PR 11)+    #[serde(default)] pub packageCooling: Vec<crate::package_cooling::Package>,     #[serde(default)] pub stackupMm: Vec<f64>,     /// each net's rail voltage for the river view; nets not named here are read from their names (GND 0, 3V3 3.3, +12V 12)     #[serde(default)] pub voltages: BTreeMap<String, f64>,@@ -48,18 +50,23 @@ pub struct Spec {     #[serde(default)] pub capCurves: Option<String>,     /// the dielectric's εr when the board's stackup does not say (4.3)     #[serde(default)] pub epsilonR: Option<f64>,+    /// the current solve's iteration budget per net (PR 5); the convergence target is unchanged+    #[serde(default = "d20000")] pub dcMaxIterations: usize,     /// set by analyze: εr per slot gap (from the board's stackup), the real gaps as given, and the real copper layers     #[serde(skip)] pub er_slots: Option<Vec<f64>>, #[serde(skip)] pub real_gaps: Vec<f64>, #[serde(skip)] pub real_copper: Vec<String>, }+fn d20000() -> usize { 20000 } fn d35() -> f64 { 35.0 } fn d60() -> f64 { 60.0 }-#[derive(Deserialize, Clone)] pub struct Load { pub amps: f64, #[serde(default)] pub maxRiseC: Option<f64> }+/// a loaded net: its current, and optionally its pins' signed currents (PR 6: positive enters the copper, negative leaves;+/// each side sums to amps), instead of the connector-to-parts heuristic+#[derive(Deserialize, Clone)] pub struct Load { pub amps: f64, #[serde(default)] pub maxRiseC: Option<f64>, #[serde(default)] pub terminals: Option<Vec<crate::terminals::Terminal>> } /// a part's dissipation: through its exposed tab's net (tabNet), or, for a leaded part with no tab, through all its leads, /// shared by pad area (issue 8) #[derive(Deserialize, Clone)] pub struct Hot { pub watts: f64, #[serde(default)] pub tabNet: String }  #[derive(Serialize, Clone)]-pub struct NetResult { pub residual: f64, pub converged: bool, pub net: String, pub amps: f64, pub layers: Vec<String>, pub peakJ: f64, pub peakAt: [f64; 2], pub peakLayer: String, pub copperMm2: f64, pub jouleW: f64, pub vias: usize, pub iterations: usize }+pub struct NetResult { pub residual: f64, pub converged: bool, pub platedPads: usize, pub platedConnections: Vec<Value>, pub net: String, pub amps: f64, pub layers: Vec<String>, pub peakJ: f64, pub peakAt: [f64; 2], pub peakLayer: String, pub copperMm2: f64, pub jouleW: f64, pub vias: usize, pub iterations: usize } #[derive(Serialize, Clone)] pub struct ChipResult { pub reference: String, pub watts: f64, pub tabNet: String, pub x: f64, pub y: f64, pub riseC: f64, pub riseMovingC: f64, pub planeW: f64, pub viaW: f64, pub vias: usize, pub otherViaW: f64, pub otherVias: usize, pub dielW: f64, pub airW: f64, pub viasWanted: usize } #[derive(Serialize, Clone)]@@ -79,6 +86,7 @@ pub struct Fields {     pub nets: Vec<NetResult>, pub chips: Vec<ChipResult>, pub issues: Vec<Issue>,     pub joule_total: f64, pub tab_total: f64,     pub heatsinks: Vec<Value>,+    pub package_cooling: Vec<Value>,     pub signals: Value,     pub heat_residual: f64, pub heat_moving_residual: f64,     /// the signed current, the solved potential and the terminals of each loaded net, for the flow and river views@@ -118,12 +126,28 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu             for r in 0..h - 1 { for c in 0..w { sys.gy[a][r * w + c] = SIGMA * t_cu[*k] * harm(f.at(r, c) as f64, f.at(r + 1, c) as f64); } }             for i in 0..n { sys.leak[a][i] = if f.v[i] > 0.0 { 1e-9 } else { 1.0 }; }         }-        let mut nvias = 0;-        for v in b.vias.iter().filter(|v| v.net == *net) {-            let Some((r, c)) = cell(v.x, v.y) else { continue };-            nvias += 1;-            let g = SIGMA * std::f64::consts::PI * ((v.drill / 2.0 + 0.025).powi(2) - (v.drill / 2.0).powi(2)) / 1.6;-            for a in 0..nl.saturating_sub(1) { sys.gz[a][r * w + c] += g; }+        // every via and plated through-hole pad barrel on this net joins its live layers, each span by its real length+        // through the stack (PR 5: a connector's plated pad joins top and bottom ground even with no via beside it)+        let nvias = b.vias.iter().filter(|v| v.net == *net).count();+        let slot = slot_gaps(spec);+        // the barrel's conductance spread around its plated wall, on cells with copper on both faces (PR 9); a wall the mesh+        // does not resolve (a small via's ring at a coarse pitch, issue 11) falls back to the nearest cells with copper on both+        // faces, so a plated hole never disconnects a net+        let barrel_cells = |v: &crate::board::Via, ka: usize, kb: usize| -> Vec<(usize, f64)> {+            let ws = crate::barrel::weights(v, b.x0, b.y0, d, w, h);+            let ok = |i: usize| fr[ka].v[i] > 0.0 && fr[kb].v[i] > 0.0;+            let norm: f64 = ws.iter().filter(|(i, _)| ok(*i)).map(|(_, q)| q).sum();+            if norm > 0.0 { return ws.into_iter().filter(|(i, _)| ok(*i)).map(|(i, q)| (i, q / norm)).collect(); }+            let Some((r0, c0)) = cell(v.x, v.y) else { return Vec::new() };+            let reach = ((v.size.max(v.drill_long.unwrap_or(0.0)) / 2.0 + d) / d).ceil() as isize;+            let mut near: Vec<(usize, f64)> = Vec::new();+            for dr in -reach..=reach { for dc in -reach..=reach { let (r, c) = (r0 as isize + dr, c0 as isize + dc); if r < 0 || c < 0 || r as usize >= h || c as usize >= w { continue; } let i = r as usize * w + c as usize; if ok(i) { near.push((i, ((dr * dr + dc * dc) as f64).sqrt())); } } }+            near.sort_by(|a, z| a.1.partial_cmp(&z.1).unwrap()); near.truncate(4); let k = near.len().max(1) as f64;+            near.into_iter().map(|(i, _)| (i, 1.0 / k)).collect()+        };+        for v in b.vias.iter().chain(b.plated_holes.iter()).filter(|v| v.net == *net) {+            let sa = SIGMA * v.barrel_area_mm2();+            for a in 0..nl.saturating_sub(1) { let (k0, k1) = (live[a], live[a + 1]); let span: f64 = slot[k0..k1].iter().sum::<f64>() + (t_cu[k0] + t_cu[k1]) / 2.0; for (i, q) in barrel_cells(v, k0, k1) { sys.gz[a][i] += sa / span.max(0.05) * q; } }         }         // terminals: the external pads (wire pads, connectors) on one side, the chips that draw the         // current on the other, each chip's share by its pad area on the net (a gate driver's sense pin@@ -150,14 +174,15 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu                 chips = if rest.is_empty() { cands.iter().cloned().filter(|p| p.reference != big).collect() } else { rest };             }         }-        if ext.is_empty() || chips.is_empty() { log(&format!("{net}: no terminals found (needs pads of two parts on the net), skipped")); continue; }+        if load.terminals.is_none() && (ext.is_empty() || chips.is_empty()) { log(&format!("{net}: no terminals found (needs pads of two parts on the net), skipped")); continue; }         let mut area_by_ref: BTreeMap<String, f64> = BTreeMap::new();         for p in &chips { *area_by_ref.entry(p.reference.clone()).or_default() += p.area.max(0.05); }         let area_total: f64 = area_by_ref.values().sum();         let mut rhs = vec![0.0; nl * n];         let mut spots_of = |p: &crate::board::Pad| -> Vec<(usize, usize, usize)> {             let Some((r0, c0)) = cell(p.x, p.y) else { return Vec::new() };-            let rad = (1.0 / d) as isize;+            // within 1 mm, or the pad's own ring when it is wider: a large plated pad's drill is a void (PR 9)+            let rad = ((1.0f64).max((p.area / std::f64::consts::PI).sqrt() + d) / d).ceil() as isize;             let mut spots: Vec<(usize, usize, usize)> = Vec::new();             for (a, k) in live.iter().enumerate() {                 if !p.layers.iter().any(|l| l == LAYERS[*k]) { continue; }@@ -170,17 +195,26 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu             spots         };         let (mut hs, mut hk) = (0usize, 0usize);-        let per_ext = load.amps / ext.len() as f64;+        // the pins and their signed currents: the spec's (PR 6), else the connectors in, the parts out by pad area+        let assignments: Vec<(&crate::board::Pad, f64)> = if let Some(rows) = &load.terminals {+            match crate::terminals::resolve(b, net, load.amps, rows) { Ok(r) => r, Err(e) => { log(&e); continue; } }+        } else {+            let mut rows: Vec<(&crate::board::Pad, f64)> = ext.iter().map(|p| (*p, load.amps / ext.len() as f64)).collect();+            for p in &chips { let share = load.amps * area_by_ref[&p.reference] / area_total / chips.iter().filter(|q| q.reference == p.reference).count() as f64; rows.push((*p, -share)); }+            rows+        };         let (mut ext_spots, mut chip_spots, mut terminals): (Vec<usize>, Vec<usize>, Vec<crate::flow::Terminal>) = (Vec::new(), Vec::new(), Vec::new());-        for p in &ext { let sp = spots_of(p); if sp.is_empty() { continue; } hs += 1; for (a, r, c) in &sp { rhs[a * n + r * w + c] += per_ext / sp.len() as f64; ext_spots.push(a * n + r * w + c); } terminals.push(crate::flow::Terminal { reference: p.reference.clone(), pad: p.pad.clone(), x: p.x, y: p.y, amps: per_ext, side: "ext".into() }); }-        let mut pinned = false; let mut pinned_at: Option<usize> = None;-        for p in &chips {-            let share = load.amps * area_by_ref[&p.reference] / area_total / chips.iter().filter(|q| q.reference == p.reference).count() as f64;-            let sp = spots_of(p); if sp.is_empty() { continue; } hk += 1;-            for (a, r, c) in &sp { rhs[a * n + r * w + c] -= share / sp.len() as f64; chip_spots.push(a * n + r * w + c); }-            terminals.push(crate::flow::Terminal { reference: p.reference.clone(), pad: p.pad.clone(), x: p.x, y: p.y, amps: share, side: "chip".into() });-            if !pinned { let (a, r, c) = sp[0]; sys.leak[a][r * w + c] = 1e3; pinned = true; pinned_at = Some(a * n + r * w + c); }+        let mut pinned = false; let mut pinned_at: Option<usize> = None; let mut missing = false;+        for (p, amps) in &assignments {+            if *amps == 0.0 { continue; }+            // an explicit pin with no copper refuses the net (PR 6); a guessed one is skipped, as before, and logged+            let sp = spots_of(p); if sp.is_empty() { if load.terminals.is_some() { missing = true; break; } log(&format!("{net}: {}.{} has no copper within 1 mm, left out of the terminals", p.reference, p.pad)); continue; }+            if *amps > 0.0 { hs += 1; } else { hk += 1; }+            for (a, r, c) in &sp { rhs[a * n + r * w + c] += amps / sp.len() as f64; if *amps > 0.0 { ext_spots.push(a * n + r * w + c); } else { chip_spots.push(a * n + r * w + c); } }+            terminals.push(crate::flow::Terminal { reference: p.reference.clone(), pad: p.pad.clone(), x: p.x, y: p.y, amps: amps.abs(), side: if *amps > 0.0 { "ext".into() } else { "chip".into() } });+            if *amps < 0.0 && !pinned { let (a, r, c) = sp[0]; sys.leak[a][r * w + c] = 1e3; pinned = true; pinned_at = Some(a * n + r * w + c); }         }+        if missing { log(&format!("{net}: a terminal has no copper under it; no partial current injection")); continue; }         if hs == 0 || hk == 0 { log(&format!("{net}: terminals not on copper (sources {hs}, sinks {hk}), skipped")); continue; }         // copper pieces with terminal current but no path to the pinned spot float on the tiny leak and stall the solve         // (issue 9: GND 6000 iterations at 1.7e-1): give each such piece its own reference@@ -195,8 +229,15 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu             if let Some(pa) = pinned_at { let root = find(&mut par, pa); let mut done: std::collections::BTreeSet<usize> = Default::default(); done.insert(root);                 for &sp in ext_spots.iter().chain(chip_spots.iter()) { let rt = find(&mut par, sp); if done.insert(rt) { let (a, i) = (sp / n, sp % n); if a < sys.nl { sys.leak[a][i] = 1e3; islands += 1; } } } }         }-        let (v, it, res) = sys.solve(&rhs, 1e-7, 20000);+        let (v, it, res) = sys.solve(&rhs, 1e-7, spec.dcMaxIterations);         log(&format!("{net}: {} A, {} layer(s), {it} iterations, residual {res:.1e}{}{}", load.amps, nl, if islands > 0 { format!(", {islands} copper piece(s) with terminals but no path to the rest, each referenced on its own") } else { String::new() }, if res > 1e-4 { " NOT CONVERGED" } else { "" }));+        // each plated hole's current between the layers it joins, and the mean density in its barrel (PR 9)+        let mut connections = Vec::new();+        for hole in b.vias.iter().chain(b.plated_holes.iter()).filter(|q| q.net == *net) {+            for a in 0..nl.saturating_sub(1) { let (k0, k1) = (live[a], live[a + 1]); let span: f64 = slot[k0..k1].iter().sum::<f64>() + (t_cu[k0] + t_cu[k1]) / 2.0; let cells = barrel_cells(hole, k0, k1); if cells.is_empty() { continue; }+                let amps = SIGMA * hole.barrel_area_mm2() / span.max(0.05) * cells.iter().map(|(i, q)| q * (v[a * n + i] - v[(a + 1) * n + i])).sum::<f64>();+                connections.push(json!({"at": [hole.x, hole.y], "from": LAYERS[k0], "to": LAYERS[k1], "amps": amps, "meanBarrelJ": amps.abs() / hole.barrel_area_mm2(), "barrelAreaMm2": hole.barrel_area_mm2(), "plated": hole.size > 0.0 && b.plated_holes.iter().any(|p| p.x == hole.x && p.y == hole.y)})); }+        }         let mut grids: Vec<Grid> = LAYERS.iter().map(|_| Grid::new(w, h)).collect();         let (mut peak, mut peak_at, mut peak_layer) = (0.0f64, [0.0, 0.0], String::new());         let mut joule = 0.0;@@ -225,7 +266,6 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu             // the signed current density and the potential on this net's copper, in the solve's own sense             // (from the connector side to the parts); flow.rs turns it into the physical sense and volts             for i in 0..n { let fi = f.v[i] as f64; if fi > 0.05 && reach[a * n + i] { let s = d * fi.max(0.25) * t_cu[*k]; nf.jx[*k].v[i] = (cx[i] / s) as f32; nf.jy[*k].v[i] = (cy[i] / s) as f32; nf.v[*k].v[i] = v[a * n + i] as f32; nf.on[*k].v[i] = 1.0; } }-            let mut vals: Vec<f32> = Vec::new();             for i in 0..n {                 let fi = f.v[i] as f64;                 if fi > 0.05 {@@ -234,14 +274,11 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu                     grids[*k].v[i] = jj as f32;                     joule += jj * jj / SIGMA * fi * d * d * t_cu[*k];                     p_joule[*k].v[i] += (jj * jj / SIGMA * fi * d * d * t_cu[*k]) as f32;-                    if fi >= 0.25 { vals.push(jj as f32); }                 }             }-            if !vals.is_empty() {-                vals.sort_by(|a, b2| a.partial_cmp(b2).unwrap());-                let cap = vals[((vals.len() as f64 * 0.995) as usize).min(vals.len() - 1)] as f64;-                for i in 0..n { if grids[*k].v[i] as f64 > cap { grids[*k].v[i] = cap as f32; } if f.v[i] >= 0.25 && (grids[*k].v[i] as f64) > peak { peak = grids[*k].v[i] as f64; peak_at = [b.x0 + (i % w) as f64 * d + d / 2.0, b.y0 + (i / w) as f64 * d + d / 2.0]; peak_layer = LAYERS[*k].to_string(); } }-            }+            // the reported peak is the modeled maximum, before the display clip (PR 8: a narrow hot spot over the limit+            // used to pass because the peak was read after the 99.5th-percentile texture clip)+            if let Some((lp, idx)) = current_peak_and_clip(&mut grids[*k], f) { if lp > peak { peak = lp; peak_at = [b.x0 + (idx % w) as f64 * d + d / 2.0, b.y0 + (idx / w) as f64 * d + d / 2.0]; peak_layer = LAYERS[*k].to_string(); } }             j_all[*k].max_with(&grids[*k]);         }         let copper: f64 = live.iter().map(|k| fr[*k].sum() * d * d).sum();@@ -252,7 +289,7 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu         let mean = |s: &[usize]| { let r: Vec<f64> = s.iter().filter(|i| reach[**i]).map(|i| v[*i]).collect(); if r.is_empty() { 0.0 } else { r.iter().sum::<f64>() / r.len() as f64 } };         nf.v_ext = mean(&ext_spots); nf.v_chip = mean(&chip_spots); nf.terminals = terminals; nf.amps = load.amps;         flow.insert(net.clone(), nf);-        nets.push(NetResult { residual: res, converged: res <= 1e-4, net: net.clone(), amps: load.amps, layers: live.iter().map(|k| LAYERS[*k].to_string()).collect(), peakJ: (peak * 10.0).round() / 10.0, peakAt: peak_at, peakLayer: peak_layer, copperMm2: (copper * 10.0).round() / 10.0, jouleW: (joule * 1000.0).round() / 1000.0, vias: nvias, iterations: it });+        nets.push(NetResult { residual: res, converged: res <= 1e-4, platedConnections: connections, platedPads: b.plated_holes.iter().filter(|v| v.net == *net).count(), net: net.clone(), amps: load.amps, layers: live.iter().map(|k| LAYERS[*k].to_string()).collect(), peakJ: peak, peakAt: peak_at, peakLayer: peak_layer, copperMm2: (copper * 10.0).round() / 10.0, jouleW: (joule * 1000.0).round() / 1000.0, vias: nvias, iterations: it });         j_net.insert(net.clone(), grids);     }     let joule_total: f64 = nets.iter().map(|x| x.jouleW).sum();@@ -289,13 +326,17 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu             for i in 0..n { sys.leak[k][i] = leak + 1e-12; }         }         for k in 0..nl - 1 { let g = K_FR4_Z * d * d / t_diel[k]; for i in 0..n { sys.gz[k][i] = g; } }-        for v in &b.vias {-            let Some((r, c)) = cell(v.x, v.y) else { continue };-            let area = std::f64::consts::PI * ((v.drill / 2.0 + 0.025).powi(2) - (v.drill / 2.0).powi(2));-            for k in 0..nl - 1 { sys.gz[k][r * w + c] += K_CU * area / t_diel[k]; }+        for v in b.vias.iter().chain(b.plated_holes.iter()) {+            let ws = crate::barrel::weights(v, b.x0, b.y0, d, w, h);+            let area = v.barrel_area_mm2();+            for k in 0..nl - 1 { for (i, q) in &ws { sys.gz[k][*i] += K_CU * area / t_diel[k].max(1e-3) * q; } }         }         let meshes:Vec<_> = spec.heatsinks.iter().map(|block| crate::heatsink::attach(&mut sys,block,b.x0,b.y0,d,h_air)).collect();-        (sys, meshes)+        let packages: Vec<_> = spec.packageCooling.iter().map(|p| {+            let nodes: Vec<_> = p.board_links.iter().filter_map(|link| { let k = LAYERS.iter().position(|l| *l == link.layer)?; let (r, c) = cell(link.at[0], link.at[1])?; Some(k * n + r * w + c) }).collect();+            p.attach(&mut sys, &nodes, h_air == H_MOVING)+        }).collect();+        (sys, meshes, packages)     };     let mut q_rhs = vec![0.0; nl * n];     for k in 0..nl { for i in 0..n {q_rhs[k*n+i] += p_joule[k].v[i] as f64;} }@@ -307,7 +348,8 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu         if cells.is_empty() { log(&format!("{reference}: {}, its watts are not placed", if leads { "no pads".to_string() } else { format!("no pad on {} (tab net)", hot.tabNet) })); continue; }         let area_sum: f64 = pads.iter().map(|p| p.area.max(1e-6)).sum();         for (pi, (x, y)) in cells.iter().enumerate() {-            let per = hot.watts * pads[pi].area.max(1e-6) / area_sum;+            // a package cooled off the board takes its watts at its junction node, not into the copper (PR 11)+            let per = if spec.packageCooling.iter().any(|p| p.reference == *reference) { 0.0 } else { hot.watts * pads[pi].area.max(1e-6) / area_sum };             let Some((r0, c0)) = cell(*x, *y) else { continue };             let rad = (1.2 / d) as isize; let mut spots = Vec::new();             for dr in -rad..=rad { for dc in -rad..=rad { let (r, c) = (r0 as isize + dr, c0 as isize + dc); if r >= 0 && c >= 0 && (r as usize) < h && (c as usize) < w && dr * dr + dc * dc <= rad * rad { spots.push((r as usize, c as usize)); } } }@@ -316,13 +358,18 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu         chip_pos.push((reference.clone(), cells[0].0, cells[0].1, hot.watts, if leads { "leads".to_string() } else { hot.tabNet.clone() }));     }     let tab_total: f64 = spec.hot.values().map(|h| h.watts).sum();-    let (sys_still, sink_meshes) = build(H_STILL);+    let (sys_still, sink_meshes, packages_still) = build(H_STILL);     q_rhs.resize(sys_still.n(),0.0);+    for (p, nodes) in spec.packageCooling.iter().zip(&packages_still) { q_rhs[nodes.junction] = spec.hot.get(&p.reference).map_or(0.0, |h| h.watts); }     let (tv, it, res) = sys_still.solve(&q_rhs, 1e-7, 8000);     log(&format!("heat, still air: {it} iterations, residual {res:.1e}"));-    let (sys_mov, _) = build(H_MOVING);+    let (sys_mov, _, packages_moving) = build(H_MOVING);     let (tm, it2, res2) = sys_mov.solve(&q_rhs, 1e-7, 8000);     log(&format!("heat, moving air: {it2} iterations"));+    let package_cooling: Vec<Value> = spec.packageCooling.iter().zip(&packages_still).zip(&packages_moving).map(|((p, a), bm)| {+        let watts = spec.hot.get(&p.reference).map_or(0.0, |h| h.watts);+        json!({"definition": p, "still": a.result(&tv, watts), "moving": bm.result(&tm, watts), "status": "supplied-parameter-steady-state-screen", "transientQualified": false, "chipRisePolicy": "maximum of board terminal region, case and junction; never hide a hot external node"})+    }).collect();     let heatsinks: Vec<Value> = spec.heatsinks.iter().zip(&sink_meshes).map(|(block,mesh)| {         let peak=|t:&[f64]|mesh.nodes.iter().map(|i|t[*i]).fold(0.0,f64::max);         let power=|t:&[f64]|mesh.contacts.iter().map(|(a,b,g)|g*(t[*a]-t[*b])).sum::<f64>();@@ -362,6 +409,8 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu         chips.push(ChipResult { reference: reference.clone(), watts: *watts, tabNet: tab.clone(), x: *x, y: *y, riseC: (tv[r0 * w + c0] * 10.0).round() / 10.0, riseMovingC: (tm[r0 * w + c0] * 10.0).round() / 10.0, planeW: (plane * 1000.0).round() / 1000.0, viaW: (via_w * 1000.0).round() / 1000.0, vias: nv, otherViaW: (other_w * 1000.0).round() / 1000.0, otherVias: no, dielW: (diel_w * 1000.0).round() / 1000.0, airW: (air_w * 1000.0).round() / 1000.0, viasWanted: via_wanted(reference) });     }     // ---- issues -------------------------------------------------------------------------------+    // a package's rise is the hottest of its board region, case and junction (PR 11)+    for (p, (a, bm)) in spec.packageCooling.iter().zip(packages_still.iter().zip(&packages_moving)) { if let Some(chip) = chips.iter_mut().find(|c| c.reference == p.reference) { chip.riseC = chip.riseC.max(tv[a.junction]).max(tv[a.case]); chip.riseMovingC = chip.riseMovingC.max(tm[bm.junction]).max(tm[bm.case]); } }     let mut issues: Vec<Issue> = Vec::new();     for x in &nets {         // IPC-2221 external, 10 C rise: about 60 A/mm2 at 35 um; a pinch well past that is a hot neck@@ -381,7 +430,7 @@ pub fn analyze(b: &Board, spec: &Spec, pitch: f64, neck: f64, log: &mut dyn FnMu         }     }     let _ = (tab_total, joule_total);-    Fields { flow, net_ids, net_names, signals:crate::signals::analyze(b,&ras,&spec.signals,&spec.copperUmByLayer,spec.copperUm), heatsinks, heat_residual:res, heat_moving_residual:res2, d, h, w, x0: b.x0, y0: b.y0, frac, j_net, cu_net, j_all, q, t, t_moving, nets, chips, issues, joule_total: (joule_total * 1000.0).round() / 1000.0, tab_total }+    Fields { package_cooling, flow, net_ids, net_names, signals:crate::signals::analyze(b,&ras,&spec.signals,&spec.copperUmByLayer,spec.copperUm), heatsinks, heat_residual:res, heat_moving_residual:res2, d, h, w, x0: b.x0, y0: b.y0, frac, j_net, cu_net, j_all, q, t, t_moving, nets, chips, issues, joule_total: (joule_total * 1000.0).round() / 1000.0, tab_total } }  /// The summary the app and the page read.@@ -399,7 +448,9 @@ pub fn summary(f: &Fields, b: &Board, spec: &Spec, textures: &[String]) -> Value     json!({         "board": {"name": spec.board, "x0": b.x0, "y0": b.y0, "widthMm": b.width(), "heightMm": b.height(), "layers": LAYERS, "copperUm": spec.copperUm, "dielectricMm":if spec.real_gaps.is_empty() {slot_gaps(spec)} else {spec.real_gaps.clone()}, "slotDielectricMm": slot_gaps(spec), "copperLayers": if spec.real_copper.is_empty() { LAYERS.iter().map(|l| l.to_string()).collect::<Vec<_>>() } else { spec.real_copper.clone() }, "copperUmByLayer":LAYERS.iter().map(|l|((*l).to_string(),spec.copperUmByLayer.get(*l).copied().unwrap_or(spec.copperUm))).collect::<BTreeMap<_,_>>(), "thicknessMm": slot_gaps(spec).iter().sum::<f64>() + (if spec.real_copper.is_empty() { LAYERS.iter().map(|l| l.to_string()).collect::<Vec<_>>() } else { spec.real_copper.clone() }).iter().map(|l|spec.copperUmByLayer.get(l).copied().unwrap_or(spec.copperUm)/1000.0).sum::<f64>(), "cellMm": f.d, "cols": f.w, "rows": f.h},         "nets": nets, "chips": chips, "issues": f.issues,-        "signals": f.signals, "heatsinks": f.heatsinks, "thermalConvergence":{"stillResidual":f.heat_residual,"movingResidual":f.heat_moving_residual},+        "peakJDefinition": {"statistic": "unclipped-cell-maximum", "minimumCopperFraction": 0.25, "displayPercentile": 0.995, "rounded": false},+        "terminalAssignments": spec.loads.iter().map(|(net, load)| (net.clone(), json!({"mode": if load.terminals.is_some() {"explicit"} else {"heuristic"}, "pins": load.terminals}))).collect::<BTreeMap<_,_>>(),+        "signals": f.signals, "heatsinks": f.heatsinks, "packageCooling": f.package_cooling, "thermalConvergence":{"stillResidual":f.heat_residual,"movingResidual":f.heat_moving_residual},         "vias": b.vias.iter().map(|v| json!({"net": v.net, "x": v.x, "y": v.y, "size": v.size, "drill": v.drill})).collect::<Vec<_>>(),         "pads": b.pads.iter().filter(|p| spec.hot.contains_key(&p.reference) || spec.loads.contains_key(&p.net)).map(|p| json!({"ref": p.reference, "pad": p.pad, "net": p.net, "x": p.x, "y": p.y})).collect::<Vec<_>>(),         "scales": {"jMax": f.j_all.iter().map(|g| g.max()).fold(0.0, f32::max), "qMax": f.q.iter().map(|g| g.max()).fold(0.0, f32::max), "tMax": tmax, "tMovingMax": f.t_moving.iter().map(|g| g.max()).fold(0.0, f32::max)},@@ -460,3 +511,43 @@ impl Clock { } #[derive(Deserialize, Clone, Debug)] pub struct Band { pub minHz: f64, pub maxHz: f64, #[serde(default)] pub ifHz: Option<f64>, #[serde(default)] pub lowSideLo: bool }++#[cfg(test)] mod terminal_tests {+    use super::*;+    // PR 6: an explicit fuse/diode path solves the same whichever way the current runs, and a sink with no copper is not silently dropped+    #[test] fn explicit_terminals() {+        let mut b = Board { x1: 10., y1: 6., layer_names: vec!["F.Cu".into(), "B.Cu".into()], copper: vec!["F.Cu".into(), "B.Cu".into()], ..Board::default() };+        b.polys.push(crate::board::Poly { net: "V".into(), layer: "F.Cu".into(), pts: vec![(0., 0.), (10., 0.), (10., 6.), (0., 6.)] });+        for (r, x) in [("F1", 1.5), ("D1", 8.5), ("R24", 5.)] { b.pads.push(crate::board::Pad { reference: r.into(), pad: "1".into(), net: "V".into(), x, y: 3., layers: vec!["F.Cu".into()], area: 1., func: String::new() }); }+        let mut spec: Spec = serde_json::from_value(json!({"copperUm": 35, "loads": {"V": {"amps": 3, "terminals": [{"reference": "F1", "pad": "1", "amps": 3}, {"reference": "D1", "pad": "1", "amps": -2.9999}, {"reference": "R24", "pad": "1", "amps": -0.0001}]}}})).unwrap();+        let f = analyze(&b, &spec, 0.5, 60., &mut |_| {});+        assert_eq!(f.nets.len(), 1); assert!(f.nets[0].converged); assert!(f.nets[0].jouleW > 0.);+        for row in spec.loads.get_mut("V").unwrap().terminals.as_mut().unwrap() { row.amps = -row.amps; }+        let rev = analyze(&b, &spec, 0.5, 60., &mut |_| {});+        assert!((f.nets[0].jouleW - rev.nets[0].jouleW).abs() < 1e-5);+        b.pads[1].x = 50.;+        assert!(analyze(&b, &spec, 0.5, 60., &mut |_| {}).nets.is_empty(), "a sink off the copper must not be silently dropped");+    }+}++/// The modeled maximum (cells at least a quarter copper) before clipping the displayed texture at its 99.5th percentile (PR 8).+fn current_peak_and_clip(grid: &mut Grid, fraction: &Grid) -> Option<(f64, usize)> {+    let mut eligible: Vec<(usize, f32)> = grid.v.iter().copied().enumerate().filter(|(i, _)| fraction.v[*i] >= 0.25).collect();+    if eligible.is_empty() { return None; }+    eligible.sort_by(|a, b| a.1.total_cmp(&b.1));+    let (index, maximum) = *eligible.last().unwrap();+    let cap = eligible[((eligible.len() as f64 * 0.995) as usize).min(eligible.len() - 1)].1;+    for v in &mut grid.v { if *v > cap { *v = cap; } }+    Some((maximum as f64, index))+}+#[cfg(test)] mod peak_tests {+    use super::*;+    #[test] fn narrow_hotspot_is_not_hidden_by_display_clip() {+        let (mut dens, mut frac) = (Grid::new(1000, 1), Grid::new(1000, 1)); dens.v.fill(20.0); frac.v.fill(1.0); dens.v[731] = 60.04;+        let (m, i) = current_peak_and_clip(&mut dens, &frac).unwrap(); assert!(m > 60.0); assert_eq!(i, 731); assert_eq!(dens.v[731], 20.0);+    }+    #[test] fn copper_fraction_cutoff() {+        let (mut dens, mut frac) = (Grid::new(2, 1), Grid::new(2, 1)); dens.v = vec![900.0, 12.0]; frac.v = vec![0.24, 0.25];+        assert_eq!(current_peak_and_clip(&mut dens, &frac), Some((12.0, 1))); frac.v.fill(0.0); assert_eq!(current_peak_and_clip(&mut dens, &frac), None);+    }+}

@@ -1,6 +1,6 @@ //! adom-fields: where the current and the heat go on your board. use sha2::{Digest, Sha256};-mod caps; mod eagle; mod emc; mod flow; mod noise; mod pdn; mod parts; mod silk; mod signals; mod heatsink; mod board; mod fields; mod raster; mod sexp; mod solve; mod tex; mod serve; mod vendor;+mod caps; mod terminals; mod barrel; mod package_cooling; mod thermal3d; mod eagle; mod emc; mod flow; mod noise; mod pdn; mod parts; mod silk; mod signals; mod heatsink; mod board; mod fields; mod raster; mod sexp; mod solve; mod tex; mod serve; mod vendor; use clap::{Parser, Subcommand}; use serde_json::{json, Value}; use std::path::PathBuf;@@ -11,6 +11,8 @@ struct Cli { #[command(subcommand)] cmd: Cmd }  #[derive(Subcommand)] enum Cmd {+    /// Solve axis-aligned 3D solids (a package, an interface, a finned sink) with explicit contacts and prescribed convection (not CFD); see docs/thermal3d.md+    Thermal3d { #[arg(long)] spec: String, #[arg(long)] out: PathBuf },     /// Solve the fields of a board: current density per loaded net, heat flow in the copper, temperature over the board; writes fields.json and the textures to --out     Analyze { #[arg(long)] board: String, #[arg(long)] spec: String, #[arg(long, default_value = "fields")] out: PathBuf, /// solve cell in mm (0.2 is a minute on a 65 x 75 mm board; 0.1 is four times longer)     #[arg(long, default_value_t = 0.2)] pitch: f64, /// the current density a loaded net may not exceed before it is an issue, A/mm² (60 is about IPC-2221 external at 35 µm for 10 °C)@@ -87,14 +89,22 @@ pub fn run_analyze(board: &str, spec_path: &str, out: &PathBuf, pitch: f64, neck     if let Some(gname) = b.pads.iter().map(|p| p.net.clone()).find(|n| crate::flow::is_ground(n)) {         if !spec.loads.contains_key(&gname) && !spec.loads.is_empty() {             let amps = spec.loads.values().map(|l| l.amps).fold(0.0, f64::max);-            if amps > 0.0 { spec.loads.insert(gname.clone(), fields::Load { amps, maxRiseC: None }); auto_ground = Some((gname, amps)); }+            if amps > 0.0 { spec.loads.insert(gname.clone(), fields::Load { amps, maxRiseC: None, terminals: None }); auto_ground = Some((gname, amps)); }         }     }     if !pitch.is_finite() || pitch <= 0.0 || pitch > b.width().min(b.height())/2.0 { return Err("pitch must be positive and resolve at least two cells across the board".into()); }     if !spec.copperUm.is_finite() || spec.copperUm<=0.0 {return Err("copperUm must be positive".into());}+    for (net, load) in &spec.loads { if let Some(rows) = &load.terminals { terminals::resolve(&b, net, load.amps, rows)?; } }+    if !(1..=100000).contains(&spec.dcMaxIterations) { return Err("dcMaxIterations must be between 1 and 100000; the convergence target is unchanged".into()); }     for (layer,t) in &spec.copperUmByLayer { if !board::LAYERS.contains(&layer.as_str()) || !t.is_finite() || *t<=0.0 {return Err("copperUmByLayer requires known copper layers and positive thicknesses".into());} }     spec.signals.validate(&spec.real_copper)?;     for hot in spec.hot.values() {if !hot.watts.is_finite() || hot.watts<0.0 {return Err("hot watts must be finite and nonnegative".into());}}+    // a package cooled off the board (PR 11): one entry per hot part, its links on that part's own pads+    for (i, package) in spec.packageCooling.iter().enumerate() {+        let cu: Vec<&str> = spec.real_copper.iter().map(String::as_str).collect();+        package.validate_board(&b, &cu, pitch)?;+        if !spec.hot.contains_key(&package.reference) || spec.packageCooling[..i].iter().any(|p| p.reference == package.reference) { return Err("packageCooling requires one entry per existing hot reference".into()); }+    }     for load in spec.loads.values() {if !load.amps.is_finite() || load.amps<0.0 {return Err("load amps must be finite and nonnegative".into());}}     if !spec.signals.capacitance.is_empty() || !spec.signals.differential_pairs.is_empty() {         for (layer_a,layer_b,distance) in spec.signals.capacitance.iter().map(|c|(&c.layer_a,&c.layer_b,c.dielectric_mm)).chain(spec.signals.differential_pairs.iter().map(|p|(&p.layer,&p.reference_layer,p.dielectric_mm))) {@@ -116,7 +126,15 @@ pub fn run_analyze(board: &str, spec_path: &str, out: &PathBuf, pitch: f64, neck     let mut lines: Vec<String> = Vec::new();     let mut f = fields::analyze(&b, &spec, pitch, neck, &mut |s| { eprintln!("{s}"); lines.push(s.to_string()); });     // the design review the caller acts on: EMC, ground, rails, the converter; every issue gets an action+    // a loaded net that was not solved: with the spec's own pins, an error, never a partial report (PR 6); with guessed pins, a+    // high issue telling the caller to name them (AI Flow refuses a solve missing a loaded net, so nothing passes on it)+    let mut unsolved: Vec<fields::Issue> = Vec::new();+    for (net, load) in &spec.loads { if load.amps > 0.0 && auto_ground.as_ref().map_or(true, |g| g.0 != *net) && !f.nets.iter().any(|r| r.net == *net) {+        if load.terminals.is_some() { return Err(format!("{net}: loaded net with explicit terminals was not solved (see the log above); no partial report is published")); }+        unsolved.push(fields::Issue { severity: "high".into(), kind: "load-unsolved".into(), what: format!("{net}: its {} A load was not solved: no terminals could be found for it (see the log)", load.amps), at: None, layer: None, net: Some(net.clone()), reference: None, fix: format!("Name {net}'s pins and their currents in spec loads.{net}.terminals (positive into the copper, negative out); a rail made inside a chip has no connector to feed it."), action: Some(serde_json::json!({"type": "review", "net": net, "why": "loaded net not solved: give spec terminals"})) });+    } }     let mut iss = std::mem::take(&mut f.issues);+    iss.extend(unsolved);     let mut review = emc::review(&f, &b, &spec, &mut iss);     let mut nz = noise::review(&f, &b, &spec, &mut iss);     let pd = pdn::review(&f, &b, &spec, &review["converter"], &mut iss, Some(out.as_path()));@@ -124,6 +142,7 @@ pub fn run_analyze(board: &str, spec_path: &str, out: &PathBuf, pitch: f64, neck     for nr in f.nets.iter().filter(|x| !x.converged) { iss.push(fields::Issue { severity: "high".into(), kind: "solve-unconverged".into(), what: format!("{}: the current solve stopped at residual {:.1e} (target 1e-7); its current density and voltage drop are not reliable", nr.net, nr.residual), at: None, layer: None, net: Some(nr.net.clone()), reference: None, fix: "Re-run at a coarser --pitch, or name this net's terminals in the spec; copper islands with terminals are referenced on their own.".into(), action: Some(json!({"type": "review", "net": nr.net, "why": "solver did not converge"})) }); }     emc::act_on_core(&mut iss, &f);     f.issues = iss;+    for package in &f.package_cooling { for state in ["still", "moving"] { let r = &package[state]; let residual = r["energyResidualW"].as_f64().unwrap_or(f64::INFINITY); let scale = 1.0 + r["sourceW"].as_f64().unwrap_or(0.0).abs() + r["boardCouplingW"].as_f64().unwrap_or(0.0).abs(); if !residual.is_finite() || residual.abs() > 1e-5 * scale { return Err(format!("package cooling energy balance failed for {} ({state})", package["definition"]["reference"])); } } }     if !f.heat_residual.is_finite() || f.heat_residual>1e-6 || !f.heat_moving_residual.is_finite() || f.heat_moving_residual>1e-6 {return Err(format!("thermal solve did not converge: still {}, moving {}; refine the model or increase the solve budget",f.heat_residual,f.heat_moving_residual));}     let tex = out.join("tex");     let mut textures: Vec<String> = Vec::new();@@ -181,13 +200,14 @@ pub fn run_analyze(board: &str, spec_path: &str, out: &PathBuf, pitch: f64, neck fn main() {     let cli = Cli::parse();     match cli.cmd {+        Cmd::Thermal3d { spec, out } => match thermal3d::run(&spec, &out) { Ok(()) => ok("3D solid thermal solve written", &[]), Err(e) => err(&e, &[]) },         Cmd::Analyze { board, spec, out, pitch, neck, no_parts, models_dir, no_silk } => {             match run_analyze(&board, &spec, &out, pitch, neck, if no_parts { None } else { Some(&models_dir) }, no_silk) {                 Ok(s) => {                     let nets = s["nets"].as_array().map(|a| a.len()).unwrap_or(0);                     let issues = s["issues"].as_array().cloned().unwrap_or_default();                     ok(&format!("{} nets solved, {} hot parts, {} issue(s), in {} s: {}", nets, s["chips"].as_array().map(|a| a.len()).unwrap_or(0), issues.len(), s["seconds"], out.join("fields.json").display()), &[-                        format!("Peak current density {:.0} A/mm²; peak rise {:.0} °C still air, {:.0} °C moving air, budget {:.0} °C; Joule loss {} W on top of {} W from the tabs.", s["scales"]["jMax"].as_f64().unwrap_or(0.0), s["heat"]["peakRiseC"].as_f64().unwrap_or(0.0), s["heat"]["peakRiseMovingC"].as_f64().unwrap_or(0.0), s["heat"]["budgetC"].as_f64().unwrap_or(0.0), s["heat"]["jouleW"], s["heat"]["tabW"]),+                        format!("Peak current density {:.0} A/mm²; peak rise {:.0} °C still air, {:.0} °C moving air, budget {:.0} °C; Joule loss {} W on top of {} W from the tabs.", s["nets"].as_array().into_iter().flatten().filter_map(|n| n["peakJ"].as_f64()).fold(0.0, f64::max), s["heat"]["peakRiseC"].as_f64().unwrap_or(0.0), s["heat"]["peakRiseMovingC"].as_f64().unwrap_or(0.0), s["heat"]["budgetC"].as_f64().unwrap_or(0.0), s["heat"]["jouleW"], s["heat"]["tabW"]),                         issues.iter().take(3).map(|i| format!("{} {}: {}", i["severity"].as_str().unwrap_or(""), i["kind"].as_str().unwrap_or(""), i["what"].as_str().unwrap_or(""))).collect::<Vec<_>>().join(" | "),                         { let r = &s["review"]; let c = &r["converter"]; let g = &r["ground"]; format!("Review: {}{}{}.", if c.is_null() { "no switching converter modelled".to_string() } else { format!("hot loop {} {} mm² (best {} mm²), emissions margin {} dB at {} MHz (FCC B, loop only), output ripple {:.1} mV", c["hotLoop"]["cap"].as_str().unwrap_or("?"), c["hotLoop"]["mm2"], c["idealMm2"], c["emissions"]["worstMarginDb"], c["emissions"]["atMHz"], c["rippleOutMv"].as_f64().unwrap_or(0.0)) }, if g.is_null() { String::new() } else { format!("; {}{} island(s), edge stitched {} %", if g["offsetMv"].is_number() { format!("ground offset {} mV, ", g["offsetMv"]) } else { "ground offset not solved (ground is not a loaded net), ".to_string() }, g["islands"].as_array().map_or(0, |x| x.len()), g["stitching"]["edgeCoveredPct"]) }, r["reference"].as_array().map(|a| a.iter().filter(|x| x["unreferencedPct"].as_f64().unwrap_or(0.0) > 5.0).take(5).map(|x| format!("; {} {}% without ground under it", x["net"].as_str().unwrap_or(""), x["unreferencedPct"])).collect::<String>()).unwrap_or_default()) },                         { let all = s["review"]["pdn"]["rails"].as_array().cloned().unwrap_or_default(); let rs: Vec<Value> = all.iter().filter(|r| r["skipped"].is_null()).cloned().collect(); let sk: String = all.iter().filter(|r| !r["skipped"].is_null()).map(|r| format!("; {} not checked: {}", r["net"].as_str().unwrap_or(""), r["skipped"].as_str().unwrap_or(""))).collect(); if rs.is_empty() { format!("PDN: no rail checked{}.", sk) } else { format!("PDN (|Z| seen from each load pin against V·tol/step): {}{}.", rs.iter().map(|r| format!("{} {}× target at {} Hz (from {})", r["net"].as_str().unwrap_or(""), r["worst"]["ratio"], r["worst"]["hz"], r["load"].as_str().unwrap_or(""))).collect::<Vec<_>>().join(", "), if rs.iter().all(|r| r["spice"]["ran"] == true) { format!("; ngspice agrees within {} %{}", rs.iter().map(|r| r["spice"]["maxErrPct"].as_f64().unwrap_or(0.0)).fold(0.0, f64::max), sk) } else { sk }) } },

@@ -0,0 +1,269 @@+//! Explicit steady-state package nodes. Heat is injected once at the junction;+//! board terminals and the external sink share it through declared resistances.+//! This is a supplied-parameter screen, not an airflow or assembly qualification.+use crate::solve::System;+use serde::{Deserialize, Serialize};+use serde_json::{json, Value};++#[derive(Clone, Deserialize, Serialize)]+#[serde(rename_all = "camelCase", deny_unknown_fields)]+pub struct Package {+    pub reference: String,+    pub junction_case_k_per_w: f64,+    pub case_sink_k_per_w: f64,+    pub sink_ambient_still_k_per_w: f64,+    pub sink_ambient_moving_k_per_w: f64,+    pub board_links: Vec<BoardLink>,+    /// Cite the part/TIM/sink and conditions, and identify assumptions.+    pub provenance: String,+}++#[derive(Clone, Deserialize, Serialize)]+#[serde(rename_all = "camelCase", deny_unknown_fields)]+pub struct BoardLink {+    pub pad: String,+    pub at: [f64; 2],+    pub layer: String,+    pub resistance_k_per_w: f64,+}++pub struct Nodes {+    pub junction: usize,+    pub case: usize,+    pub sink: usize,+    pub links: Vec<(usize, f64)>,+    pub sink_conductance: f64,+}++impl Package {+    pub fn validate_board(+        &self,+        board: &crate::board::Board,+        layers: &[&str],+        pitch: f64,+    ) -> Result<(), String> {+        self.validate()?;+        if !pitch.is_finite() || pitch <= 0.0 {+            return Err("invalid package grid pitch".into());+        }+        let width = ((board.width() / (pitch / 2.0)).round() as usize).max(1) / 2;+        let height = ((board.height() / (pitch / 2.0)).round() as usize).max(1) / 2;+        for link in &self.board_links {+            if !layers.contains(&link.layer.as_str()) {+                return Err(format!(+                    "{}: package link names an unsupported layer",+                    self.reference+                ));+            }+            let matches: Vec<_> = board+                .pads+                .iter()+                .filter(|p| {+                    p.reference == self.reference+                        && p.pad == link.pad+                        && (p.x - link.at[0]).abs() < 1e-6+                        && (p.y - link.at[1]).abs() < 1e-6+                        && p.layers.contains(&link.layer)+                })+                .collect();+            if matches.len() != 1 {+                return Err(format!(+                    "{}: package link must identify exactly one physical pad on its copper layer",+                    self.reference+                ));+            }+            if link.at[0] < board.x0+                || link.at[0] >= board.x0 + width as f64 * pitch+                || link.at[1] < board.y0+                || link.at[1] >= board.y0 + height as f64 * pitch+            {+                return Err(format!(+                    "{}: package link lies outside thermal grid",+                    self.reference+                ));+            }+        }+        Ok(())+    }++    pub fn validate(&self) -> Result<(), String> {+        if self.reference.trim().is_empty() || self.provenance.trim().is_empty() {+            return Err("packageCooling requires reference and parameter provenance".into());+        }+        for r in [+            self.junction_case_k_per_w,+            self.case_sink_k_per_w,+            self.sink_ambient_still_k_per_w,+            self.sink_ambient_moving_k_per_w,+        ]+        .into_iter()+        .chain(self.board_links.iter().map(|l| l.resistance_k_per_w))+        {+            if !r.is_finite() || r <= 0.0 || !(1.0 / r).is_finite() {+                return Err(format!(+                    "{}: package resistances must be finite and positive",+                    self.reference+                ));+            }+        }+        for (i, link) in self.board_links.iter().enumerate() {+            if link.pad.trim().is_empty() || !link.at.iter().all(|x| x.is_finite()) {+                return Err(format!(+                    "{}: every board link needs a pad and finite coordinates",+                    self.reference+                ));+            }+            if self.board_links[..i]+                .iter()+                .any(|p| p.pad == link.pad && p.at == link.at && p.layer == link.layer)+            {+                return Err(format!("{}: duplicate package board link", self.reference));+            }+        }+        Ok(())+    }++    /// Caller resolves and validates physical pads and layer indices first.+    pub fn attach(&self, sys: &mut System, board_nodes: &[usize], moving: bool) -> Nodes {+        assert_eq!(board_nodes.len(), self.board_links.len());+        let board_count = sys.nl * sys.h * sys.w;+        assert!(board_nodes.iter().all(|i| *i < board_count));+        let sink_conductance = 1.0+            / if moving {+                self.sink_ambient_moving_k_per_w+            } else {+                self.sink_ambient_still_k_per_w+            };+        let junction = sys.add_node(0.0);+        let case = sys.add_node(0.0);+        let sink = sys.add_node(sink_conductance);+        sys.edges+            .push((junction, case, 1.0 / self.junction_case_k_per_w));+        sys.edges.push((case, sink, 1.0 / self.case_sink_k_per_w));+        let links = board_nodes+            .iter()+            .zip(&self.board_links)+            .map(|(node, link)| {+                let conductance = 1.0 / link.resistance_k_per_w;+                sys.edges.push((case, *node, conductance));+                (*node, conductance)+            })+            .collect();+        Nodes {+            junction,+            case,+            sink,+            links,+            sink_conductance,+        }+    }+}++impl Nodes {+    pub fn result(&self, temperatures: &[f64], watts: f64) -> Value {+        let board_w: f64 = self+            .links+            .iter()+            .map(|(i, g)| g * (temperatures[self.case] - temperatures[*i]))+            .sum();+        let ambient_w = temperatures[self.sink] * self.sink_conductance;+        json!({+            "junctionRiseC": temperatures[self.junction], "caseRiseC": temperatures[self.case],+            "sinkRiseC": temperatures[self.sink], "boardCouplingW": board_w,+            "sinkToAmbientW": ambient_w, "sourceW": watts,+            "energyResidualW": watts - board_w - ambient_w,+            "boardLinks": self.links.iter().map(|(i,g)| json!({+                "boardRiseC": temperatures[*i],+                "caseToBoardW": g*(temperatures[self.case]-temperatures[*i])+            })).collect::<Vec<_>>(),+        })+    }+}++#[cfg(test)]+mod tests {+    use super::*;+    fn package() -> Package {+        Package {+            reference: "Q1".into(),+            junction_case_k_per_w: 0.5,+            case_sink_k_per_w: 1.0,+            sink_ambient_still_k_per_w: 3.0,+            sink_ambient_moving_k_per_w: 2.0,+            board_links: vec![],+            provenance: "Analytical test resistances, not physical qualification".into(),+        }+    }+    fn close(a: f64, b: f64) {+        assert!((a - b).abs() < 1e-8, "{a} != {b}");+    }+    #[test]+    fn isolated_series_path_has_exact_temperatures_and_conserves_power() {+        let p = package();+        let mut s = System::new(1, 1, 1);+        s.leak[0][0] = 1.0;+        let n = p.attach(&mut s, &[], false);+        let mut q = vec![0.0; s.n()];+        q[n.junction] = 2.0;+        let (t, _, res) = s.solve(&q, 1e-12, 100);+        assert!(res < 1e-10);+        close(t[n.sink], 6.0);+        close(t[n.case], 8.0);+        close(t[n.junction], 9.0);+        close(n.result(&t, 2.0)["energyResidualW"].as_f64().unwrap(), 0.0);+    }+    #[test]+    fn parallel_board_path_splits_power_without_double_injection() {+        let mut p = package();+        p.board_links.push(BoardLink {+            pad: "2".into(),+            at: [0.0, 0.0],+            layer: "F.Cu".into(),+            resistance_k_per_w: 2.0,+        });+        let mut s = System::new(1, 1, 1);+        s.leak[0][0] = 0.5;+        let n = p.attach(&mut s, &[0], false);+        let mut q = vec![0.0; s.n()];+        q[n.junction] = 2.0;+        let (t, _, _) = s.solve(&q, 1e-12, 100);+        let r = n.result(&t, 2.0);+        close(t[n.case], 4.0);+        close(t[n.junction], 5.0);+        close(t[0], 2.0);+        close(r["boardCouplingW"].as_f64().unwrap(), 1.0);+        close(r["sinkToAmbientW"].as_f64().unwrap(), 1.0);+    }+    #[test]+    fn hotter_board_can_send_heat_back_into_external_sink() {+        let mut p = package();+        p.board_links.push(BoardLink {+            pad: "2".into(),+            at: [0.0, 0.0],+            layer: "F.Cu".into(),+            resistance_k_per_w: 2.0,+        });+        let mut s = System::new(1, 1, 1);+        s.leak[0][0] = 0.5;+        let n = p.attach(&mut s, &[0], false);+        let mut q = vec![0.0; s.n()];+        q[0] = 10.0;+        q[n.junction] = 1.0;+        let (t, _, _) = s.solve(&q, 1e-12, 100);+        let r = n.result(&t, 1.0);+        assert!(r["boardCouplingW"].as_f64().unwrap() < 0.0);+        close(r["energyResidualW"].as_f64().unwrap(), 0.0);+        close(t[0] * 0.5 + r["sinkToAmbientW"].as_f64().unwrap(), 11.0);+    }+    #[test]+    fn invalid_or_untraceable_resistances_are_rejected() {+        for r in [0.0, -1.0, f64::NAN, f64::INFINITY] {+            let mut p = package();+            p.case_sink_k_per_w = r;+            assert!(p.validate().is_err());+        }+        let mut p = package();+        p.provenance.clear();+        assert!(p.validate().is_err());+    }+}

@@ -0,0 +1,98 @@+//! Explicit signed pin currents. Positive enters copper; negative leaves it.+use crate::board::{Board, Pad};+use serde::{Deserialize, Serialize};+use std::collections::BTreeSet;++#[derive(Clone, Debug, Deserialize, Serialize)]+#[serde(deny_unknown_fields)]+pub struct Terminal {+    pub reference: String,+    pub pad: String,+    pub amps: f64,+    /// Optional board-space pad center in mm, for repeated logical pin numbers.+    #[serde(default, skip_serializing_if = "Option::is_none")]+    pub at: Option<[f64; 2]>,+}++pub fn resolve<'a>(board: &'a Board, net: &str, amps: f64, rows: &[Terminal])+    -> Result<Vec<(&'a Pad, f64)>, String>+{+    if rows.is_empty() { return Err(format!("{net}: explicit terminals must not be empty")); }+    let mut seen = BTreeSet::new();+    let mut result = Vec::new();+    let (mut positive, mut negative) = (0.0, 0.0);+    for t in rows {+        if !t.amps.is_finite() || t.at.is_some_and(|xy| xy.iter().any(|v| !v.is_finite())) {+            return Err(format!("{net}: terminal currents must be finite and pin assignments unique"));+        }+        let pads: Vec<_> = board.pads.iter().enumerate().filter(|(_, p)| {+            p.reference == t.reference && p.pad == t.pad+                && t.at.is_none_or(|[x, y]| (p.x - x).abs() <= 1e-6 && (p.y - y).abs() <= 1e-6)+        }).collect();+        if pads.len() != 1 || pads[0].1.net != net {+            return Err(format!("{net}: {}.{} must identify exactly one pad on this net", t.reference, t.pad));+        }+        if !seen.insert(pads[0].0) {+            return Err(format!("{net}: physical pad assigned more than once"));+        }+        positive += t.amps.max(0.0);+        negative += (-t.amps).max(0.0);+        result.push((pads[0].1, t.amps));+    }+    let tolerance = 1e-9_f64.max(amps.abs() * 1e-8);+    if !amps.is_finite() || amps < 0.0 || (positive - amps).abs() > tolerance || (negative - amps).abs() > tolerance {+        return Err(format!("{net}: sources ({positive} A) and sinks ({negative} A) must each equal load ({amps} A)"));+    }+    Ok(result)+}++#[cfg(test)]+mod tests {+    use super::*;+    fn board() -> Board {+        Board { pads: [("J1", "1", "V"), ("F1", "1", "V"), ("R1", "1", "V"), ("J2", "1", "GND")].into_iter().map(|(r,p,n)| Pad { reference:r.into(), pad:p.into(), net:n.into(), x:0., y:0., layers:vec!["F.Cu".into()], area:1., func:String::new() }).collect(), ..Board::default() }+    }+    fn row(reference: &str, amps: f64) -> Terminal { Terminal { reference:reference.into(), pad:"1".into(), amps, at: None } }+    #[test] fn main_and_sense_branches_preserve_assigned_currents() {+        let b=board();+        let r=resolve(&b,"V",3.,&[row("J1",3.),row("F1",-2.9999),row("R1",-0.0001)]).unwrap();+        assert_eq!(r[1].0.reference,"F1"); assert_eq!(r[2].1,-0.0001);+    }+    #[test] fn reject_wrong_net_missing_duplicate_and_unbalanced_terminals() {+        let b=board();+        for rows in [vec![],vec![row("J1",3.),row("J2",-3.)],vec![row("J1",3.),row("MISSING",-3.)],vec![row("J1",3.),row("F1",-2.)],vec![row("J1",3.),row("F1",-1.5),row("F1",-1.5)],vec![row("J1",f64::NAN),row("F1",-3.)]] {+            assert!(resolve(&b,"V",3.,&rows).is_err());+        }+    }+    #[test] fn reject_ambiguous_repeated_pad_number() {+        let mut b=board(); b.pads.push(b.pads[0].clone());+        assert!(resolve(&b,"V",3.,&[row("J1",3.),row("F1",-3.)]).is_err());+    }+    #[test] fn repeated_pin_physical_pads_require_explicit_centers() {+        let mut b=board();+        let mut second=b.pads[1].clone(); second.x=2.; b.pads.push(second);+        let mut first=row("F1",-1.5); first.at=Some([0.,0.]);+        let mut second=row("F1",-1.5); second.at=Some([2.,0.]);+        let resolved=resolve(&b,"V",3.,&[row("J1",3.),first.clone(),second.clone()]).unwrap();+        assert_eq!(resolved[1].0.x,0.); assert_eq!(resolved[2].0.x,2.);+        assert!(resolve(&b,"V",3.,&[row("J1",3.),first.clone(),first.clone()]).is_err());+        for at in [Some([1.,0.]),Some([f64::NAN,0.]),None] {+            second.at=at;+            assert!(resolve(&b,"V",3.,&[row("J1",3.),first.clone(),second.clone()]).is_err());+        }+    }+    #[test] fn selector_cannot_bypass_wrong_net_or_duplicate_assignment() {+        let b=board(); let mut wrong=row("J2",-3.);wrong.at=Some([0.,0.]);+        assert!(resolve(&b,"V",3.,&[row("J1",3.),wrong]).is_err());+        let mut selected=row("F1",-1.5);selected.at=Some([0.,0.]);+        assert!(resolve(&b,"V",3.,&[row("J1",3.),row("F1",-1.5),selected]).is_err());+    }+    #[test] fn old_json_remains_valid_and_new_selector_roundtrips() {+        let old:Terminal=serde_json::from_str(r#"{"reference":"F1","pad":"1","amps":-3}"#).unwrap();+        assert!(old.at.is_none());+        let mut new=old;new.at=Some([111.138,90.08]);+        let restored:Terminal=serde_json::from_str(&serde_json::to_string(&new).unwrap()).unwrap();+        assert_eq!(restored.at,new.at);+    }++}\ No newline at end of file

@@ -0,0 +1,434 @@+//! Steady 3D solid conduction. Prescribed convection is NOT a CFD fan model.+use crate::solve::System;+use serde::{Deserialize, Serialize};+use serde_json::{json, Value};+use sha2::{Digest, Sha256};+use std::{+    collections::{HashMap, HashSet},+    path::Path,+};++#[derive(Clone, Deserialize, Serialize)]+#[serde(rename_all = "camelCase", deny_unknown_fields)]+pub struct Body {+    pub name: String,+    /// Union of axis-aligned boxes [xmin,ymin,zmin,xmax,ymax,zmax], mm.+    pub boxes_mm: Vec<[f64; 6]>,+    pub conductivity_w_mk: f64,+    /// Exposed-face h, ordered -x,+x,-y,+y,-z,+z. Zero means insulated.+    pub convection_w_m2k: [f64; 6],+    pub watts: f64,+    pub provenance: String,+}+#[derive(Clone, Deserialize, Serialize)]+#[serde(rename_all = "camelCase", deny_unknown_fields)]+pub struct Contact {+    pub bodies: [String; 2],+    pub resistance_m2k_w: f64,+    pub provenance: String,+}+#[derive(Clone, Deserialize, Serialize)]+#[serde(rename_all = "camelCase", deny_unknown_fields)]+pub struct Spec {+    pub name: String,+    pub ambient_c: f64,+    pub pitch_mm: f64,+    pub bodies: Vec<Body>,+    pub contacts: Vec<Contact>,+    pub assumptions: Vec<String>,+}+#[derive(Clone, Serialize)]+#[serde(rename_all = "camelCase")]+struct Cell {+    body: usize,+    bounds_mm: [f64; 6],+    temperature_c: f64,+}+fn positive(x: f64) -> bool {+    x.is_finite() && x > 0.0+}+fn nonnegative(x: f64) -> bool {+    x.is_finite() && x >= 0.0+}+fn axis(mut edges: Vec<f64>, pitch: f64) -> Result<Vec<f64>, String> {+    edges.sort_by(f64::total_cmp);+    edges.dedup_by(|a, b| (*a - *b).abs() < 1e-8);+    let mut out = Vec::new();+    for p in edges.windows(2) {+        let n = ((p[1] - p[0]) / pitch).ceil() as usize;+        if n > 10000 || out.len() + n > 10000 {+            return Err("axis exceeds 10000 cells".into());+        }+        for i in 0..n {+            out.push(p[0] + (p[1] - p[0]) * i as f64 / n as f64);+        }+    }+    out.push(*edges.last().ok_or("empty geometry")?);+    Ok(out)+}+fn near_index(a: &[f64], x: f64) -> usize {+    a.iter().position(|v| (*v - x).abs() < 1e-7).unwrap()+}++pub fn solve(spec: &Spec) -> Result<Value, String> {+    if !spec.ambient_c.is_finite()+        || spec.ambient_c <= -273.15+        || !positive(spec.pitch_mm)+        || spec.bodies.is_empty()+        || spec.assumptions.is_empty()+    {+        return Err(+            "require ambient above absolute zero, positive pitch, bodies and assumptions".into(),+        );+    }+    let mut names = HashMap::new();+    let mut axes: [Vec<f64>; 3] = Default::default();+    for (i, b) in spec.bodies.iter().enumerate() {+        if b.name.trim().is_empty()+            || names.insert(b.name.clone(), i).is_some()+            || b.provenance.trim().is_empty()+            || b.boxes_mm.is_empty()+            || !positive(b.conductivity_w_mk)+            || !nonnegative(b.watts)+            || b.convection_w_m2k.iter().any(|h| !nonnegative(*h))+        {+            return Err("invalid or duplicate body; require material/geometry provenance".into());+        }+        for bb in &b.boxes_mm {+            for d in 0..3 {+                if !bb[d].is_finite() || !bb[d + 3].is_finite() || bb[d + 3] - bb[d] < 1e-7 {+                    return Err("invalid box bounds".into());+                }+                axes[d].extend([bb[d], bb[d + 3]]);+            }+        }+    }+    let mut contacts = HashMap::new();+    for c in &spec.contacts {+        let a = *names.get(&c.bodies[0]).ok_or("unknown contact body")?;+        let b = *names.get(&c.bodies[1]).ok_or("unknown contact body")?;+        if a == b+            || !nonnegative(c.resistance_m2k_w)+            || c.provenance.trim().is_empty()+            || contacts+                .insert((a.min(b), a.max(b)), c.resistance_m2k_w)+                .is_some()+        {+            return Err("invalid/duplicate contact".into());+        }+    }+    for a in &mut axes {+        *a = axis(std::mem::take(a), spec.pitch_mm)?;+    }+    let dims = [axes[0].len() - 1, axes[1].len() - 1, axes[2].len() - 1];+    let total = dims+        .iter()+        .try_fold(1usize, |a, b| a.checked_mul(*b))+        .ok_or("grid overflow")?;+    if total > 2_000_000 {+        return Err("grid exceeds 2 million cells; simplify geometry or increase pitch".into());+    }+    let flat = |p: [usize; 3]| p[0] + dims[0] * (p[1] + dims[1] * p[2]);+    let mut owner = vec![usize::MAX; total];+    for (bi, b) in spec.bodies.iter().enumerate() {+        for bb in &b.boxes_mm {+            let lo = [+                near_index(&axes[0], bb[0]),+                near_index(&axes[1], bb[1]),+                near_index(&axes[2], bb[2]),+            ];+            let hi = [+                near_index(&axes[0], bb[3]),+                near_index(&axes[1], bb[4]),+                near_index(&axes[2], bb[5]),+            ];+            for z in lo[2]..hi[2] {+                for y in lo[1]..hi[1] {+                    for x in lo[0]..hi[0] {+                        let o = &mut owner[flat([x, y, z])];+                        if *o != usize::MAX && *o != bi {+                            return Err(format!(+                                "overlapping bodies {} and {}",+                                spec.bodies[*o].name, b.name+                            ));+                        }+                        *o = bi;+                    }+                }+            }+        }+    }+    let mut map = vec![usize::MAX; total];+    let mut cells = Vec::new();+    let mut coords = Vec::new();+    let mut volumes = vec![0.0; spec.bodies.len()];+    let mut sys = System::new(0, 1, 1);+    for z in 0..dims[2] {+        for y in 0..dims[1] {+            for x in 0..dims[0] {+                let p = [x, y, z];+                let o = owner[flat(p)];+                if o == usize::MAX {+                    continue;+                }+                let bb = [+                    axes[0][x],+                    axes[1][y],+                    axes[2][z],+                    axes[0][x + 1],+                    axes[1][y + 1],+                    axes[2][z + 1],+                ];+                volumes[o] += (bb[3] - bb[0]) * (bb[4] - bb[1]) * (bb[5] - bb[2]);+                map[flat(p)] = sys.add_node(0.0);+                coords.push(p);+                cells.push(Cell {+                    body: o,+                    bounds_mm: bb,+                    temperature_c: 0.0,+                });+            }+        }+    }+    let mut used = HashSet::new();+    let mut exposed = Vec::new();+    let mut edge_meta = Vec::new();+    for (i, c) in cells.iter().enumerate() {+        let p = coords[i];+        let size = [+            c.bounds_mm[3] - c.bounds_mm[0],+            c.bounds_mm[4] - c.bounds_mm[1],+            c.bounds_mm[5] - c.bounds_mm[2],+        ];+        for d in 0..3 {+            for sign in 0..2 {+                let mut nb = p;+                let j = if (sign == 0 && p[d] == 0) || (sign == 1 && p[d] + 1 == dims[d]) {+                    usize::MAX+                } else {+                    if sign == 0 {+                        nb[d] -= 1;+                    } else {+                        nb[d] += 1;+                    }+                    map[flat(nb)]+                };+                let area = size[(d + 1) % 3] * size[(d + 2) % 3] * 1e-6;+                if j == usize::MAX {+                    let h = spec.bodies[c.body].convection_w_m2k[2 * d + sign];+                    // Includes half-cell conduction to the convecting face.+                    let g = if h > 0.0 {+                        area / (1.0 / h + size[d] * 0.0005 / spec.bodies[c.body].conductivity_w_mk)+                    } else {+                        0.0+                    };+                    sys.extra_leak[i] += g;+                    exposed.push((i, d, sign, area, g));+                } else if sign == 1 {+                    let other = &cells[j];+                    let key = (c.body.min(other.body), c.body.max(other.body));+                    let rc = if c.body == other.body {+                        0.0+                    } else {+                        used.insert(key);+                        *contacts.get(&key).ok_or_else(|| {+                            format!(+                                "touching bodies require explicit contact: {} / {}",+                                spec.bodies[c.body].name, spec.bodies[other.body].name+                            )+                        })?+                    };+                    let g = area+                        / (size[d] * 0.0005 / spec.bodies[c.body].conductivity_w_mk+                            + (other.bounds_mm[d + 3] - other.bounds_mm[d]) * 0.0005+                                / spec.bodies[other.body].conductivity_w_mk+                            + rc);+                    sys.edges.push((i, j, g));+                    edge_meta.push((d, area));+                }+            }+        }+    }+    if contacts.keys().any(|k| !used.contains(k)) {+        return Err("declared contact has no shared face; check placement".into());+    }+    // Every connected solid island must have a path to ambient, even unpowered islands.+    let mut adj = vec![Vec::new(); cells.len()];+    for &(a, b, _) in &sys.edges {+        adj[a].push(b);+        adj[b].push(a);+    }+    let mut seen = vec![false; cells.len()];+    for root in 0..cells.len() {+        if !seen[root] {+            let mut stack = vec![root];+            seen[root] = true;+            let mut sink = false;+            while let Some(a) = stack.pop() {+                sink |= sys.extra_leak[a] > 0.0;+                for &b in &adj[a] {+                    if !seen[b] {+                        seen[b] = true;+                        stack.push(b);+                    }+                }+            }+            if !sink {+                return Err("solid island has no thermal path to ambient".into());+            }+        }+    }+    let rhs: Vec<f64> = cells+        .iter()+        .map(|c| {+            let b = c.bounds_mm;+            spec.bodies[c.body].watts * (b[3] - b[0]) * (b[4] - b[1]) * (b[5] - b[2])+                / volumes[c.body]+        })+        .collect();+    let (rise, iterations, residual) = sys.solve(&rhs, 1e-10, 20000);+    let input: f64 = rhs.iter().sum();+    let rejected: f64 = rise.iter().zip(&sys.extra_leak).map(|(t, g)| t * g).sum();+    let balance = (input - rejected).abs() / input.max(1e-12);+    if !residual.is_finite()+        || residual > 1e-8+        || rise.iter().any(|v| !v.is_finite() || *v < -1e-6)+        || balance > 1e-6+    {+        return Err(format!(+            "failed convergence/energy balance: residual {residual}, balance {balance}"+        ));+    }+    for (c, t) in cells.iter_mut().zip(&rise) {+        c.temperature_c = spec.ambient_c + t;+    }+    let body_results:Vec<Value>=spec.bodies.iter().enumerate().map(|(i,b)|{+        let cs:Vec<&Cell>=cells.iter().filter(|c|c.body==i).collect();+        json!({"name":b.name,"watts":b.watts,"minC":cs.iter().map(|c|c.temperature_c).fold(f64::INFINITY,f64::min),"maxC":cs.iter().map(|c|c.temperature_c).fold(f64::NEG_INFINITY,f64::max),"volumeMm3":volumes[i]})+    }).collect();+    let faces:Vec<Value>=exposed.iter().map(|&(i,d,s,area,g)|json!({"cell":i,"axis":d,"positive":s==1,"areaM2":area,"wattsToAmbient":g*rise[i]})).collect();+    let flux:Vec<Value>=sys.edges.iter().zip(&edge_meta).map(|(&(a,b,g),&(d,area))|json!({"from":a,"to":b,"axis":d,"watts":g*(rise[a]-rise[b]),"fluxWm2":g*(rise[a]-rise[b])/area})).collect();+    Ok(+        json!({"schema":"adom-fields-thermal3d-v1","name":spec.name,"method":"Steady finite-volume solid conduction; prescribed convection; no CFD, radiation or heat-pipe transport","ambientC":spec.ambient_c,"pitchMm":spec.pitch_mm,"assumptions":spec.assumptions,"input":spec,"gridDimensions":dims,"cellCount":cells.len(),"iterations":iterations,"relativeResidual":residual,"inputWatts":input,"rejectedWatts":rejected,"relativeEnergyError":balance,"bodies":body_results,"cells":cells,"exposedFaces":faces,"conductiveFaces":flux}),+    )+}++pub fn run(path: &str, out: &Path) -> Result<(), String> {+    let bytes = std::fs::read(path).map_err(|e| e.to_string())?;+    let spec: Spec = serde_json::from_slice(&bytes).map_err(|e| e.to_string())?;+    let mut result = solve(&spec)?;+    result["inputSha256"] = json!(format!("{:x}", Sha256::digest(&bytes)));+    result["solverVersion"] = json!(env!("CARGO_PKG_VERSION"));+    if std::fs::read(path).map_err(|e| e.to_string())? != bytes {+        return Err("input changed during solve".into());+    }+    std::fs::create_dir_all(out).map_err(|e| e.to_string())?;+    std::fs::write(+        out.join("thermal3d.json"),+        serde_json::to_vec(&result).unwrap(),+    )+    .map_err(|e| e.to_string())?;+    Ok(())+}++#[cfg(test)]+mod tests {+    use super::*;+    fn slab() -> Spec {+        Spec {+            name: "slab".into(),+            ambient_c: 20.0,+            pitch_mm: 1.0,+            bodies: vec![Body {+                name: "solid".into(),+                boxes_mm: vec![[0., 0., 0., 10., 10., 10.]],+                conductivity_w_mk: 10.,+                convection_w_m2k: [0., 0., 0., 0., 0., 100.],+                watts: 1.,+                provenance: "analytical fixture".into(),+            }],+            contacts: vec![],+            assumptions: vec!["uniform volumetric heat, only top convects".into()],+        }+    }+    #[test]+    fn slab_matches_analytic_and_balance() {+        let s = slab();+        let r = solve(&s).unwrap();+        let max = r["bodies"][0]["maxC"].as_f64().unwrap(); // T(bottom cell center) = Ta+Q/hA+q'''*(L²-z²)/2k+        let expected = 20. + 100. + 1000000. * (0.01f64.powi(2) - 0.0005f64.powi(2)) / 20.;+        assert!((max - expected).abs() < 0.03);+        assert!(r["relativeEnergyError"].as_f64().unwrap() < 1e-8);+    }+    #[test]+    fn insulated_island_rejected() {+        let mut s = slab();+        s.bodies[0].convection_w_m2k = [0.; 6];+        assert!(solve(&s).unwrap_err().contains("island"));+    }+    #[test]+    fn overlap_rejected() {+        let mut s = slab();+        let mut b = s.bodies[0].clone();+        b.name = "other".into();+        s.bodies.push(b);+        assert!(solve(&s).unwrap_err().contains("overlapping"));+    }+    #[test]+    fn contact_drop_and_missing_contact() {+        let mut s = slab();+        s.bodies[0].convection_w_m2k = [0.; 6];+        let mut b = s.bodies[0].clone();+        b.name = "sink".into();+        b.boxes_mm = vec![[0., 0., 10., 10., 10., 20.]];+        b.watts = 0.;+        b.convection_w_m2k[5] = 100.;+        s.bodies.push(b);+        assert!(solve(&s).unwrap_err().contains("explicit contact"));+        s.contacts.push(Contact {+            bodies: ["solid".into(), "sink".into()],+            resistance_m2k_w: 0.,+            provenance: "fixture".into(),+        });+        let a = solve(&s).unwrap();+        s.contacts[0].resistance_m2k_w = 0.001;+        let b = solve(&s).unwrap();+        assert!(+            (b["bodies"][0]["maxC"].as_f64().unwrap()+                - a["bodies"][0]["maxC"].as_f64().unwrap()+                - 10.)+                .abs()+                < 1e-5+        );+    }+    #[test]+    fn wrong_contact_rejected() {+        let mut s = slab();+        s.contacts.push(Contact {+            bodies: ["solid".into(), "missing".into()],+            resistance_m2k_w: 0.,+            provenance: "test".into(),+        });+        assert!(solve(&s).is_err());+    }+    #[test]+    fn zero_power_is_ambient() {+        let mut s = slab();+        s.bodies[0].watts = 0.;+        let r = solve(&s).unwrap();+        assert_eq!(r["bodies"][0]["maxC"], 20.);+    }+    #[test]+    fn refinement_converges() {+        let mut s = slab();+        s.pitch_mm = 2.;+        let a = solve(&s).unwrap();+        s.pitch_mm = 0.5;+        let b = solve(&s).unwrap();+        assert!(+            (a["bodies"][0]["maxC"].as_f64().unwrap() - b["bodies"][0]["maxC"].as_f64().unwrap())+                .abs()+                < 0.1+        );+    }+}