Handboek · hoofdstuk 36Numerieke methoden en wetenschappelijk rekenen

Een NumPy-array heeft een vaste dtype, vorm (shape) en geheugenindeling. t=np.arange(0,3601,60) geeft 61 tijdstippen in seconde. T=5+20*np.exp(-t/10800) evalueert dezelfde formule over alle punten. NumPy voert veel elementbewerkingen in gecompileerde code uit; de uitdrukking maakt ook zichtbaar dat alle samples dezelfde regel volgen.

Een tabel met n tijden en twee temperaturen heeft vorm (n,2). De vector met twee capaciteiten heeft vorm (2,); broadcasting deelt elke kolom door de bijbehorende capaciteit. Vorm (n,1) naast (n,) kan onbedoeld een (n, n)-array maken. Controleer vormen, eenheden en assen voordat een snelle berekening een groot verkeerd resultaat produceert. Gebruik @ voor matrixproducten, * voor elementgewijze producten.

Een nulpunt voldoet aan $f(x)=0$. Als f continu is en $f(a)f(b)<0$, ligt minstens één nulpunt tussen a en b. Bereken het midden en behoud de helft met tekenwisseling. De invariant is dus een continu interval met een ingesloten wortel. Bij een sprong kan tekenwisseling zonder nulpunt bestaan; continuïteit is essentieel.

Na k stappen is de intervalbreedte $(b-a)/2^k$. Geef je het midden terug, dan is de fout hoogstens $(b-a)/2^{k+1}$. Dit levert een echte garantie zonder afgeleide. Stop zowel op een absolute als relatieve breedte in praktische code; een kleine functiewaarde alleen is onvoldoende bij een vlakke functie.

Lineariseer bij $x_k$: $f(x)\approx f(x_k)+f^{\prime }(x_k)(x-x_k)$. Zet dit gelijk aan nul:

Bij een eenvoudige wortel, gladde f en een goede start is de convergentie lokaal kwadratisch: het aantal goede cijfers kan ongeveer verdubbelen. De methode heeft geen algemene globale garantie. Een kleine afgeleide geeft een enorme sprong; een slecht startpunt kan naar een andere wortel leiden of blijven slingeren. Voor $f=x^3-2x+2$ gaat start 0 naar 1 en weer naar 0.

Een safeguarded methode bewaart een betrouwbaar bracket en aanvaardt Newton alleen als de stap daarbinnen ligt en voldoende voortgang maakt; anders volgt bisectie. In meerdere dimensies wordt de afgeleide een Jacobiaan en los je een lineair stelsel op. Dat vraagt goed geschaalde variabelen en een controle op singulariteit.

Het model voorspelt ${\hat{y}}_i(\theta )$; de data geven $y_i$. Kleinste kwadraten kiest parameters die $J(\theta )={\sum }_i[({\hat{y}}_i-y_i)/{\sigma }_i{]}^2$ minimaliseren. Gewichten passen bij de meetonzekerheid. Ze maken een grof instrument niet belangrijker doordat het veel samples levert. Gecorreleerde samples vragen een covariantiemodel of passende subsampling.

Een warmteproef moet input en beginvoorwaarden kennen. De proef hier gebruikt de modulecapaciteiten $C_a=80kJ/K$, $C_b=1,10MJ/K$, buiten 5 °C, koppeling 30 W/K en verlies 20,2 W/K. De capaciteiten worden vastgehouden; de fit schat alleen verlies UA en koppeling H. Twee sensoren en opwarm- én afkoelfase maken die parameters beter onderscheidbaar. Dit model hoort bij de module; de 540 kJ/K van e39 is een afzonderlijke snelle stationsknoop.

De lucht ontvangt vermogen P, verliest warmte naar buiten en wisselt warmte uit met de buffer. De buffer ontvangt precies die laatste stroom. Daarom:

Optellen schrapt de interne uitwisseling: $C_a{\dot{T}}_a+C_b{\dot{T}}_b=P-UA(T_a-T_o)$. Dat is een controle op energiebewaring. Een fout in het teken van H zou bij twee gelijke temperaturen onzichtbaar kunnen blijven maar later energie creëren. Integratie en residuen worden met dezelfde modelroutine berekend, zodat grafiek en fit niet uiteenlopen.

Euler gebruikt de helling aan het begin: $y_{k+1}=y_k+hf(y_k,t_k)$. De lokale fout per stap is $O(h^2)$, de globale fout over een vaste tijd $O(h)$. RK4 gebruikt hellingen aan begin, twee middenpunten en einde: $k_1=f(y,t)$, $k_2=f(y+h{k}_1/2,t+h/2)$, $k_3=f(y+h{k}_2/2,t+h/2)$, $k_4=f(y+h{k}_3,t+h)$. De combinatie is:

De gewichten vallen niet uit een gemiddelde-intuïtie alleen af te leiden; ze zijn gekozen om de Taylorontwikkeling tot en met h⁴ te laten overeenkomen. De globale fout is bij voldoende gladheid $O(h^4)$. Bij een inputschakelaar is het proces niet glad over de stap: laat de integrator precies op het schakelmoment eindigen. Adaptieve methoden schatten de fout en verkleinen de stap waar nodig.

Voor $\dot{y}=\lambda y$ geeft Euler $y_{k+1}=(1+h\lambda )y_k$. Een werkelijk uitdovend systeem blijft alleen numeriek uitdoven als $|1+h\lambda |<1$. Bij $\lambda =-1/\tau$ betekent dit $0<h<2\tau$. Als $h>\tau$ wisselt de numerieke oplossing van teken, hoewel echte afkoeling dat niet doet. Voor monotone afkoeling is h≤τ noodzakelijk in deze test; voor nauwkeurigheid kies je veel kleiner.

In het tweeknoopsmodel bepalen beide eigenwaarden van de systeemmatrix de stapgrens. De snelle lucht kan een stap klein maken terwijl de buffer urenlang verandert: een stijve situatie. RK4 heeft een groter maar eindig stabiliteitsgebied. Een impliciete methode kan helpen, maar verplaatst het werk naar stelseloplossingen. Onvoorwaardelijke stabiliteit betekent evenmin automatisch nauwkeurigheid.

Float 64 bewaart ongeveer 16 significante decimale cijfers. Ronding ontstaat bij elke niet exact representeerbare operatie. $0,1+0,2$ hoeft binair niet exact 0,3 te zijn. Gebruik een tolerantie die uit de toepassing volgt; willekeurig ‘meer decimalen’ tonen is geen onzekerheidsanalyse.

Bij subtractie van bijna gelijke getallen verdwijnt relatieve informatie. $1-\cos x$ is voor klein x nauwkeuriger als $2{\sin }^2(x/2)$. Voor exponentiële kleine veranderingen helpt np.expm 1(x). Numeriek differentiëren $(y(t+h)-y(t))/h$ verlaagt truncatiefout bij kleiner h, maar versterkt meetruis en afronding door 1/h. Een goede stap is dus een afweging. Conditionering betreft het probleem; stabiliteit betreft het algoritme dat het oplost.

Een klein residu kan ontstaan door te veel vrije parameters. Fit niet alle warmtecapaciteiten, koppelingen en sensorvertragingen tegelijk op één korte luchtmeting. Plot residuen tegen tijd, temperatuur en input; een langdurig tekenpatroon wijst op vertraging of een ontbrekende knoop. Valideer op een andere vermogenscyclus en buitenconditie.

Bereken gevoeligheden numeriek: verander één parameter licht en bekijk de verandering van de voorspelling. Bij bijna parallelle gevoeligheidskolommen zijn parameters onderling uitwisselbaar. Een covariance uit lokale linearizatie is dan instabiel. Rapporteer de vastgehouden parameters, grenzen, meetruis, validatiefout en geldigheidsgebied. Een ‘fit’ op synthetische data controleert het programma; ze toont nog niet dat Habitat werkelijk dit model volgt.