Energiberegning etter NS 3031:2025 med timesmetode og energirapport

enok leser en bygningsbeskrivelse i TOML, henter typisk klimaår for kommunen,
simulerer 8760 timer etter NS-EN ISO 52016-1 og skriver energirapport, CSV og JSON.

Motor
- Femnodemodell for opake bygningsdeler, tonodemodell for vinduer, luftnode og
  strålingsnode per sone, implisitt løsning per time med egen tett LU-løser.
- Varme- og kjølebehov ved eksakt interpolasjon mellom to løsninger, siden systemet
  er lineært. Vinduslufting for bygninger uten kjøleanlegg.
- Solberegning etter ISO 52010-1, kontrollert mot pvlib med avvik under 1 prosent.

Beregningsgang
- Beregningspunktene A til E: netto, brutto, tilført, levert per energibærer og
  vektet levert energi, med solceller, klimakorrigering og effektbehov.
- Energikarakter etter Enovas skala fra 1.1.2026, med arealkorreksjon.

Inndata
- TOML med scenarier for tiltak, og import fra IFC, gbXML og regneark.
- Normerte verdier, karakterskala, vektingsfaktorer og klimadata i data/,
  som kan overstyres med --data.

Kontroll
- Stasjonært varmetap mot U·A·ΔT innenfor 2 prosent, energibalansen innenfor
  2 prosent, varmetapsbudsjettet for TEK17-referansebygget identisk med publisert
  tabell. Kjent avvik i romoppvarming mot Simien 8 er dokumentert i README.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_015GeH3y8fwTaTudaD88P57o

src/motor/iso52016.rs

@@ -0,0 +1,518 @@
+//! Dynamisk timesberegning av netto energibehov (beregningspunkt A) etter
+//! NS-EN ISO 52016-1, slik NS 3031:2025 krever.
+//!
+//! Hver opake bygningsdel modelleres med fem noder der varmekapasiteten fordeles etter
+//! masseplassering (ISO 52016-1 tabell B.13/B.14), og vinduer med to noder uten masse.
+//! Sonen har én luftnode med varmekapasitet for luft og inventar, og én strålingsnode
+//! («stjernenode») som fordeler langbølget strålingsutveksling mellom innvendige overflater.
+//! Systemet løses implisitt (bakover-Euler) for hver time.
+//!
+//! Varme- og kjølebehovet bestemmes som i standarden: systemet løses først uten tilført
+//! effekt og deretter med en prøveeffekt. Siden systemet er lineært, gir interpolasjon
+//! mellom de to løsningene nøyaktig den effekten som treffer settpunktet.
+
+use crate::klima::{sol, Vaerdata};
+use crate::modell::lost::{Lost, Regulering};
+use crate::modell::Mot;
+use crate::motor::{lin, ventilasjon, RHO_C_LUFT};
+use crate::{Result, TIMER};
+use serde::Serialize;
+
+/// Konvektiv varmeovergang innvendig overflate, W/(m²K). Gir sammen med h_ri R_si = 0,13.
+const H_CI: f64 = 2.5;
+/// Langbølget strålingsutveksling innvendig, W/(m²K).
+const H_RI: f64 = 5.13;
+/// Konvektiv varmeovergang utvendig, W/(m²K).
+const H_CE: f64 = 20.0;
+/// Langbølget stråling utvendig, W/(m²K). Gir sammen med h_ce R_se = 0,04.
+const H_RE: f64 = 4.14;
+/// Konvektiv andel av interne varmetilskudd (ISO 52016-1).
+const KONVEKTIV_ANDEL: f64 = 0.4;
+/// Innstråling der solskjerming trer i kraft, W/m².
+const SKJERMING_TERSKEL: f64 = 300.0;
+/// Utetemperatur der manuell solskjerming antas tatt i bruk, °C.
+const SKJERMING_TEMP: f64 = 10.0;
+const DT: f64 = 3600.0;
+/// Prøveeffekt per m² ved bestemmelse av varme-/kjølebehov.
+const PROVEEFFEKT: f64 = 10.0;
+
+/// Nodeoppsett for én bygningsdel.
+struct Del {
+    areal: f64,
+    /// Indeks til delens første node i tilstandsvektoren.
+    forste: usize,
+    /// Antall noder (5 for opake deler, 2 for vinduer og dører).
+    antall: usize,
+    /// Konduktans mellom nabonoder, W/K (allerede multiplisert med areal).
+    h_indre: f64,
+    /// Konduktans mot ytre miljø fra node 1, W/K.
+    h_ytre: f64,
+    /// Varmekapasitet per node, J/K.
+    kapasitet: Vec<f64>,
+    glass: bool,
+    mot: Mot,
+    b_faktor: f64,
+    absorptans: f64,
+    helning: f64,
+    /// Solinnstråling på flaten, W/m², for hele året.
+    innstraaling: Vec<f64>,
+    g_effektiv: f64,
+    skjermingsfaktor: f64,
+    automatisk_skjerming: bool,
+    horisontfaktor: f64,
+}
+
+/// Resultatet av en årssimulering: alle timesserier, W (energistrømmer) og °C.
+#[derive(Serialize, Clone, Debug)]
+pub struct Timeserie {
+    pub t_ute: Vec<f64>,
+    pub t_inne: Vec<f64>,
+    pub t_operativ: Vec<f64>,
+    /// Romoppvarming, post 1a, W.
+    pub romoppvarming: Vec<f64>,
+    /// Ventilasjonsvarme, post 1b, W.
+    pub ventilasjonsvarme: Vec<f64>,
+    /// Romkjøling, post 3a, W.
+    pub romkjoling: Vec<f64>,
+    /// Ventilasjonskjøling, post 3b, W.
+    pub ventilasjonskjoling: Vec<f64>,
+    /// Varmtvann, post 2, W.
+    pub varmtvann: Vec<f64>,
+    /// Vifter, post 4a, W.
+    pub vifter: Vec<f64>,
+    /// Pumper, post 4b, W.
+    pub pumper: Vec<f64>,
+    /// Belysning, post 5, W.
+    pub belysning: Vec<f64>,
+    /// Teknisk utstyr, post 6, W.
+    pub utstyr: Vec<f64>,
+    /// Sum solinnstråling gjennom vinduer, W.
+    pub soltilskudd: Vec<f64>,
+    /// Interne varmetilskudd, W.
+    pub internlast: Vec<f64>,
+    /// Solstråling på horisontal flate, W/m².
+    pub sol_horisontal: Vec<f64>,
+    /// Luftmengde ved vinduslufting, m³/h.
+    pub vinduslufting: Vec<f64>,
+}
+
+impl Timeserie {
+    fn ny() -> Self {
+        let n = || vec![0.0; TIMER];
+        Timeserie {
+            t_ute: n(),
+            t_inne: n(),
+            t_operativ: n(),
+            romoppvarming: n(),
+            ventilasjonsvarme: n(),
+            romkjoling: n(),
+            ventilasjonskjoling: n(),
+            varmtvann: n(),
+            vifter: n(),
+            pumper: n(),
+            belysning: n(),
+            utstyr: n(),
+            soltilskudd: n(),
+            internlast: n(),
+            sol_horisontal: n(),
+            vinduslufting: n(),
+        }
+    }
+
+    /// Årssum av en serie, kWh.
+    pub fn sum_kwh(serie: &[f64]) -> f64 {
+        serie.iter().sum::<f64>() / 1000.0
+    }
+
+    /// Månedssummer, kWh.
+    pub fn maaned_kwh(serie: &[f64]) -> [f64; 12] {
+        let mut ut = [0.0; 12];
+        for (t, v) in serie.iter().enumerate() {
+            ut[crate::maaned_for_time(t)] += v / 1000.0;
+        }
+        ut
+    }
+}
+
+/// Uregulerbart varmetilskudd i sonen fra systemtap, W.
+///
+/// Tap fra tappevannsbereder og varmtvannsrør i oppvarmet sone varmer rommene og reduserer
+/// romoppvarmingsbehovet, samtidig som energien telles under tappevann. Distribusjonstap for
+/// romvarme i oppvarmet sone regnes ikke med her: den varmen ville uansett blitt avgitt i sonen,
+/// og påvirker derfor bare fordelingen mellom rør og varmegiver, ikke energibehovet.
+#[derive(Clone, Copy, Debug, Default)]
+pub struct Systemtilskudd {
+    /// Konstant tilskudd fra tappevannsanlegget, W.
+    pub tappevann: f64,
+}
+
+pub fn simuler(l: &Lost, v: &Vaerdata, tilskudd: Systemtilskudd) -> Result<Timeserie> {
+    let deler = bygg_deler(l, v);
+    let n_noder: usize = 2 + deler.iter().map(|d| d.antall).sum::<usize>();
+    let areal_innvendig: f64 = deler.iter().map(|d| d.areal).sum();
+    let c_luft = l.innvendig_varmekapasitet * l.bra;
+    let h_inf = ventilasjon::infiltrasjonsluftmengde(l) * RHO_C_LUFT;
+    let bakketemp = bakketemperatur(v);
+
+    let mut ts = Timeserie::ny();
+    let mut tilstand = vec![l.settpunkt_oppvarming; n_noder];
+    let mut t_a_forrige = l.settpunkt_oppvarming;
+
+    // Innkjøring: gjenta de første to ukene til starttilstanden har stabilisert seg.
+    for innkjoring in 0..3 {
+        let skriv = innkjoring == 2;
+        let start_tilstand = tilstand.clone();
+        let timer: Box<dyn Iterator<Item = usize>> = if skriv {
+            Box::new(0..TIMER)
+        } else {
+            Box::new(0..(24 * 14))
+        };
+        if !skriv {
+            tilstand = start_tilstand;
+        }
+        for t in timer {
+            let klokke = t % 24;
+            let t_ute = v.t_ute[t];
+            let t_himmel = sol::himmeltemperatur(v, t);
+
+            // ── laster ──────────────────────────────────────────────
+            let bel = l.laster.belysning[klokke] * l.bra;
+            let utstyr_el = l.laster.utstyr_effekt[klokke] * l.bra;
+            let utstyr_varme = l.laster.utstyr_varme[klokke] * l.bra;
+            let personer = l.laster.personer[klokke] * l.bra;
+            let vv = l.laster.varmtvann[klokke] * l.bra;
+            let pumper = l.laster.pumper[klokke] * l.bra;
+            let vifter = ventilasjon::vifteeffekt(l, klokke);
+            // Vifteeffekten ender som varme i lufta (ISO 52016-1: hele SFP-effekten).
+            let internlast = bel + utstyr_varme + personer + tilskudd.tappevann + vifter;
+
+            // ── sol ─────────────────────────────────────────────────
+            let mut sol_gjennom = 0.0;
+            for d in &deler {
+                if d.glass && d.g_effektiv > 0.0 {
+                    sol_gjennom += soltransmisjon(d, t, t_ute);
+                }
+            }
+
+            // ── ventilasjon ─────────────────────────────────────────
+            let q_vent = ventilasjon::luftmengde(l, klokke);
+            let h_ve = q_vent * RHO_C_LUFT;
+            let t_sup_hx = ventilasjon::tilluftstemperatur(l, t_ute, t_a_forrige);
+            let (t_sup, vent_varme) =
+                if l.ventilasjon.varmebatteri && t_sup_hx < l.ventilasjon.tilluft_settpunkt {
+                    (
+                        l.ventilasjon.tilluft_settpunkt,
+                        h_ve * (l.ventilasjon.tilluft_settpunkt - t_sup_hx),
+                    )
+                } else {
+                    (t_sup_hx, 0.0)
+                };
+
+            // ── settpunkt ───────────────────────────────────────────
+            let sett_varme = if l.laster.oppvarming_driftstid.aktiv(klokke) {
+                l.settpunkt_oppvarming
+            } else {
+                l.settpunkt_natt
+            };
+
+            // ── løsning ─────────────────────────────────────────────
+            let prove = PROVEEFFEKT * l.bra;
+            let ctrl = |s: &[f64]| -> f64 {
+                match l.regulering {
+                    Regulering::Luft => s[0],
+                    Regulering::Operativ => 0.5 * s[0] + 0.5 * s[1],
+                }
+            };
+            let lag_ctx = |h_luft: f64| Kontekst {
+                deler: &deler,
+                n_noder,
+                areal_innvendig,
+                c_luft,
+                h_inf,
+                h_ve: h_ve + h_luft,
+                // Lufteluften kommer rett utenfra, uten varmegjenvinning.
+                t_sup: if h_ve + h_luft > 1e-9 {
+                    (h_ve * t_sup + h_luft * t_ute) / (h_ve + h_luft)
+                } else {
+                    t_sup
+                },
+                t_ute,
+                t_himmel,
+                t_bakke: bakketemp[crate::maaned_for_time(t)],
+                t_a_forrige,
+                internlast,
+                sol_gjennom,
+                time: t,
+            };
+
+            // Vinduslufting: øk luftmengden til den operative temperaturen holder settpunktet.
+            let mut h_luft = 0.0;
+            let mut ctx = lag_ctx(0.0);
+            let mut los0 = ctx.los(&tilstand, 0.0)?;
+            if let Some((sett_luft, maks_luft)) = l.lufting {
+                if ctrl(&los0) > sett_luft {
+                    let h_maks = maks_luft * l.bra * RHO_C_LUFT;
+                    let ctx_maks = lag_ctx(h_maks);
+                    let los_maks = ctx_maks.los(&tilstand, 0.0)?;
+                    if ctrl(&los_maks) >= sett_luft {
+                        h_luft = h_maks;
+                        ctx = ctx_maks;
+                        los0 = los_maks;
+                    } else {
+                        let (mut lav, mut hoy) = (0.0, h_maks);
+                        for _ in 0..12 {
+                            let midt = 0.5 * (lav + hoy);
+                            let k = lag_ctx(midt);
+                            let l0 = k.los(&tilstand, 0.0)?;
+                            if ctrl(&l0) > sett_luft {
+                                lav = midt;
+                            } else {
+                                hoy = midt;
+                            }
+                            h_luft = midt;
+                            ctx = k;
+                            los0 = l0;
+                        }
+                    }
+                }
+            }
+            let los1 = ctx.los(&tilstand, prove)?;
+            let c0 = ctrl(&los0);
+            let c1 = ctrl(&los1);
+            let stigning = (c1 - c0) / prove; // K per W
+
+            let mut fi = 0.0;
+            if c0 < sett_varme - 1e-9 && stigning > 1e-12 {
+                fi = (sett_varme - c0) / stigning;
+                if let Some(maks) = l.maks_varmeeffekt {
+                    fi = fi.min(maks);
+                }
+            } else if let Some(sett_kjol) = l.settpunkt_kjoling {
+                if c0 > sett_kjol + 1e-9 && stigning > 1e-12 {
+                    fi = (sett_kjol - c0) / stigning;
+                    if let Some(maks) = l.maks_kjoleeffekt {
+                        fi = fi.max(-maks);
+                    }
+                }
+            }
+            // Eksakt tilstand ved interpolasjon (systemet er lineært).
+            let andel = fi / prove;
+            for i in 0..n_noder {
+                tilstand[i] = los0[i] + andel * (los1[i] - los0[i]);
+            }
+            t_a_forrige = tilstand[0];
+
+            if skriv {
+                ts.t_ute[t] = t_ute;
+                ts.t_inne[t] = tilstand[0];
+                ts.t_operativ[t] = 0.5 * tilstand[0] + 0.5 * tilstand[1];
+                ts.romoppvarming[t] = fi.max(0.0);
+                ts.romkjoling[t] = (-fi).max(0.0);
+                ts.ventilasjonsvarme[t] = vent_varme;
+                ts.varmtvann[t] = vv;
+                ts.vifter[t] = vifter;
+                ts.pumper[t] = pumper;
+                ts.belysning[t] = bel;
+                ts.utstyr[t] = utstyr_el;
+                ts.soltilskudd[t] = sol_gjennom;
+                ts.internlast[t] = internlast;
+                ts.sol_horisontal[t] = v.ghi[t];
+                ts.vinduslufting[t] = h_luft / RHO_C_LUFT;
+            }
+        }
+    }
+    Ok(ts)
+}
+
+struct Kontekst<'a> {
+    deler: &'a [Del],
+    n_noder: usize,
+    areal_innvendig: f64,
+    c_luft: f64,
+    h_inf: f64,
+    h_ve: f64,
+    t_sup: f64,
+    t_ute: f64,
+    t_himmel: f64,
+    t_bakke: f64,
+    t_a_forrige: f64,
+    internlast: f64,
+    sol_gjennom: f64,
+    time: usize,
+}
+
+impl Kontekst<'_> {
+    /// Setter opp og løser timesligningene med gitt tilført effekt til sonen.
+    fn los(&self, forrige: &[f64], fi_hc: f64) -> Result<Vec<f64>> {
+        let n = self.n_noder;
+        let mut a = vec![0.0; n * n];
+        let mut b = vec![0.0; n];
+        let idx = |i: usize, j: usize| i * n + j;
+
+        // ── rad 0: sonens luftnode ──────────────────────────────────
+        let c_dt = self.c_luft / DT;
+        a[idx(0, 0)] = c_dt + self.h_ve + self.h_inf;
+        b[0] = c_dt * forrige[0]
+            + self.h_ve * self.t_sup
+            + self.h_inf * self.t_ute
+            + KONVEKTIV_ANDEL * self.internlast
+            + fi_hc;
+        for d in self.deler {
+            let indre = d.forste + d.antall - 1;
+            let h = d.areal * H_CI;
+            a[idx(0, 0)] += h;
+            a[idx(0, indre)] -= h;
+        }
+
+        // ── rad 1: strålingsnoden (ingen masse, netto null) ─────────
+        for d in self.deler {
+            let indre = d.forste + d.antall - 1;
+            let h = d.areal * H_RI;
+            a[idx(1, 1)] += h;
+            a[idx(1, indre)] -= h;
+        }
+        b[1] = 0.0;
+
+        // ── bygningsdelene ──────────────────────────────────────────
+        let rad_last = (1.0 - KONVEKTIV_ANDEL) * self.internlast + self.sol_gjennom;
+        for d in self.deler {
+            let t_ytre = self.ytre_temperatur(d);
+            let sol_absorbert = if d.glass {
+                0.0
+            } else {
+                d.areal * d.absorptans * d.innstraaling[self.time]
+            };
+            for k in 0..d.antall {
+                let i = d.forste + k;
+                let c = d.kapasitet[k] / DT;
+                a[idx(i, i)] += c;
+                b[i] += c * forrige[i];
+                if k == 0 {
+                    a[idx(i, i)] += d.h_ytre;
+                    b[i] += d.h_ytre * t_ytre + sol_absorbert;
+                }
+                if k + 1 < d.antall {
+                    a[idx(i, i)] += d.h_indre;
+                    a[idx(i, i + 1)] -= d.h_indre;
+                }
+                if k > 0 {
+                    a[idx(i, i)] += d.h_indre;
+                    a[idx(i, i - 1)] -= d.h_indre;
+                }
+                if k == d.antall - 1 {
+                    a[idx(i, i)] += d.areal * (H_CI + H_RI);
+                    a[idx(i, 0)] -= d.areal * H_CI;
+                    a[idx(i, 1)] -= d.areal * H_RI;
+                    b[i] += rad_last * d.areal / self.areal_innvendig;
+                }
+            }
+        }
+        lin::los(&mut a, &mut b, n)?;
+        Ok(b)
+    }
+
+    fn ytre_temperatur(&self, d: &Del) -> f64 {
+        match d.mot {
+            Mot::Grunn => self.t_bakke,
+            Mot::Uoppvarmet => self.t_a_forrige - d.b_faktor * (self.t_a_forrige - self.t_ute),
+            Mot::Ute => {
+                // Himmelandel etter flatens helning; resten ser bakken med utelufttemperatur.
+                let f_himmel = (1.0 + d.helning.to_radians().cos()) / 2.0;
+                let t_str = f_himmel * self.t_himmel + (1.0 - f_himmel) * self.t_ute;
+                (H_CE * self.t_ute + H_RE * t_str) / (H_CE + H_RE)
+            }
+        }
+    }
+}
+
+fn soltransmisjon(d: &Del, t: usize, t_ute: f64) -> f64 {
+    let i_sol = d.innstraaling[t] * d.horisontfaktor;
+    let skjerm = if d.skjermingsfaktor < 1.0
+        && i_sol > SKJERMING_TERSKEL
+        && (d.automatisk_skjerming || t_ute > SKJERMING_TEMP)
+    {
+        d.skjermingsfaktor
+    } else {
+        1.0
+    };
+    d.areal * d.g_effektiv * skjerm * i_sol
+}
+
+fn bygg_deler(l: &Lost, v: &Vaerdata) -> Vec<Del> {
+    let mut deler = Vec::with_capacity(l.flater.len());
+    let mut neste = 2usize;
+    for f in &l.flater {
+        let antall = if f.glass { 2 } else { 5 };
+        let r_total = 1.0 / f.u_verdi.max(1e-6);
+        let r_konstruksjon = (r_total - 1.0 / (H_CI + H_RI) - 1.0 / (H_CE + H_RE)).max(1e-3);
+        let gap = (antall - 1) as f64;
+        let h_indre = f.areal * gap / r_konstruksjon;
+        let h_ytre = f.areal * (H_CE + H_RE);
+        let kapasitet = fordel_kapasitet(f, antall);
+        let innstraaling = if f.glass || f.mot == Mot::Ute {
+            sol::innstraaling_aar(v, f.helning, f.orientering, l.albedo)
+        } else {
+            vec![0.0; TIMER]
+        };
+        deler.push(Del {
+            areal: f.areal,
+            forste: neste,
+            antall,
+            h_indre,
+            h_ytre,
+            kapasitet,
+            glass: f.glass,
+            mot: f.mot,
+            b_faktor: f.b_faktor,
+            absorptans: f.absorptans,
+            helning: f.helning,
+            innstraaling,
+            g_effektiv: f.g_effektiv,
+            skjermingsfaktor: f.skjermingsfaktor,
+            automatisk_skjerming: f.automatisk_skjerming,
+            horisontfaktor: (1.0 - f.horisont / 90.0).clamp(0.0, 1.0),
+        });
+        neste += antall;
+    }
+    deler
+}
+
+/// Fordeler varmekapasiteten over nodene etter ISO 52016-1 tabell B.13.
+fn fordel_kapasitet(f: &crate::modell::lost::LostFlate, antall: usize) -> Vec<f64> {
+    use crate::modell::Masseplassering::*;
+    if antall == 2 {
+        return vec![0.0, 0.0];
+    }
+    let k = f.kappa * f.areal;
+    let mut c = vec![0.0; 5];
+    match f.masseplassering {
+        Innvendig => c[4] = k,
+        Utvendig => c[0] = k,
+        Begge => {
+            c[0] = k / 2.0;
+            c[4] = k / 2.0;
+        }
+        Midt => c[2] = k,
+        Fordelt => {
+            c[0] = k / 8.0;
+            c[1] = k / 4.0;
+            c[2] = k / 4.0;
+            c[3] = k / 4.0;
+            c[4] = k / 8.0;
+        }
+    }
+    c
+}
+
+/// Dempet bakketemperatur per måned etter NS-EN ISO 13370 (forenklet): årsmiddel med
+/// halvert amplitude og én måneds forsinkelse.
+fn bakketemperatur(v: &Vaerdata) -> [f64; 12] {
+    let aar = v.aarsmiddeltemperatur();
+    let maaned = v.maanedsmiddeltemperatur();
+    let mut ut = [0.0; 12];
+    for m in 0..12 {
+        let forsinket = maaned[(m + 11) % 12];
+        ut[m] = aar + 0.5 * (forsinket - aar);
+    }
+    ut
+}