![](https://www.ntbinfo.no/data/images/00304/f6275cea-f278-4251-a592-ff9f54ba0352.jpg)
Koblede diffligninger
IN-KJM1900, femte forelesning, 13.november
![](https://www.ntbinfo.no/data/images/00304/f6275cea-f278-4251-a592-ff9f54ba0352.jpg)
Mål for økta:¶
- Litt repetisjon om interpolasjon
- Fortsettelse p? deloppgave 2
- Hydrologisk delmodell
- Kjemisk delmodell
- Ny teori:
- Koblede differensialligninger
Orakeltjeneste, uke 46
Orakelhjelp tilbys ved Kjemisk Institutt p? f?lgende tidspunkter denne uka
Tirsdag 6/11 | 14.00-16.00 |
Onsdag 7/11 | 13.00-15.00 |
Torsdag 8/11 | 13.00-15.00 |
Innlevering av deloppgave 2
- Fristen for ? levere deloppgave 2 av prosjektet er s?ndag 18. november klokka 23.59.
- Oppgaven leveres (som normalt) p? Devilry.? ?
- Du oppfordres sterkt til ? benytte gruppetimer og samretting p? torsdag klokka 10.15-12.00
- Sp?rsm?l om prosjektet?
Interpolasjon (repetisjon)
Interpolasjon er en teknikk hvor vi benytter diskrete datapunkter for ? gj?re et mer eller mindre kvalifisert gjett p? hva som finnes mellom. Dette er p? en m?te en "anti-diskretisering".Interpolasjon
I dette kurset vil vi benytte en interpolasjons-funksjon fra scipy, som importeres og benyttes p? f?lgende m?te:from scipy.interpolate import interp1d
import numpy as np
# data punkter for testing
data_x = np.array([1,2,3,4])
data_y = np.array([1,7,3,5])
y_interp = interp1d(data_x, data_y, 3) #gj?r en kubisk (3. ordens) interpolasjon
x = np.linspace(1,4,100) #?k oppl?sningen p? grid
y = y_interp(x) #beregn funksjonsverdiene til den interpolerte funksjonen
Eksamensrelevante sp?rsm?l
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/bmod.png)
Prosjektets del 2: Birkenesmodellen
- Numerisk modell for vannbevegelse og kjemi i jordv?ske.
- Tilpasset observasjoner i nedb?rsfelt i Birkenes kommune, Aust-Agder.
- Flere publikasjoner p? 70-, 80- og 90-tallet.
- Kjemisk modell koblet til hydrologisk modell.
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/birkenes-revidert.png)
Hydrologisk delmodell: detaljer
- Nedb?rsfelt modellert som reservoarer $\mathbf{A}$ og $\mathbf{B}$
- For et reservoar, for eksempel merket $\mathbf{M}$:
- Vannstand som funksjon av tid: $M(t)$
- Str?m ut: $Q_M(t)$
- Evaporasjon: $E_M(t)$
- Nedb?r ("precipitation") $P_M(t)$
- Generelt: subindeks indikerer reservoar.
- Avrenning til bekk. Fysiske ma?linger gjort her.
-
Her benyttes fet font for "merkelappene" p? reservoarene.
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/at.png)
Reservoar $\mathbf{A}$
Vi beskriver f?rst modellen for reservoar A: Vannstand som funksjon av tid: $A(t)$ m?les i $mm$, hvor $mm$ tilsvarer vannmengde per flateenhet slik at $1mm = 1 \frac{l}{m^2}$. Vann str?mmer inn i reservoaret kun gjennom nedb?ren $P(t)$. Vann str?mmer ut fra reservoaret gjennom evaporasjon $E_A$ og en bekk $Q_A$. Vannstr?m/fluks m?les i $mm/dag = mm \cdot dag^{-1}$.![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/at.png)
Reservoar $\mathbf{A}$
Regler for reservoar $\mathbf{A}$:- Ingen avrenning ($Q_A = 0$) n?r vannstanden $A(t)$ er mindre enn $A_{min}$
- Avrenning intensiveres n?r vannstanden $A(t)$ ?ker over $A_{min}$
- Evaporasjon er en funksjon av temperatur $E_A := E_A(T)$
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/at.png)
Reservoar $\mathbf{A}$
Vannstanden i $\mathbf{A}$ p?virkes positivt (?kes) av nedb?ren $P$, og negativt av evapotranspirasjon $E$ og utstr?mming til bekken $Q_A$, alts? er $A := A(t, P_A, E_A, Q_A)$.- Symbolet $\delta$ er en liten gresk delta ($\Delta$), og brukes konvensjonelt for ? beskrive endringer og differanser. Andre steder kan vi bruke $\Delta$ eller noen ganger $h$.
- Dimensjonen til $\delta t$ er tid, vannstr?m/fluks har dimensjon vannmengde per flateenhet per tid ($mm/d?gn$) og $\delta A$ har dimensjonen vannmengde per flateenhet.
- Sp?rsm?l til diskusjon: g?r dimensjonene i uttrykket for endringen opp? Diskuter.
$$ mm = \frac{mm}{d?gn} \cdot \text{d?gn} + \frac{mm}{d?gn} \cdot \text{d?gn} - \frac{mm}{d?gn} \cdot \text{d?gn} $$
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/at.png)
Str?m $Q_A$ som funksjon av vannstanden $A(t)$
N?r vannstanden $A(t)$ er h?y er str?mmen ut av $\mathbf{A}$ sterk, og n?r vannstanden er lav er str?mmen svak. Str?mmen er derfor en funksjon av vannstanden $A(t)$, og i modellen v?r uttrykkes dette som![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/at.png)
Evapotranspirasjon $E_A$ som funksjon av temperaturen $T$
Ved h?ye temperaturer g?r mange prosesser fortere. I dette tilfellet har vi h?yere fordamping og ?kt opptak i vegetasjon, og vi uttrykker det som en line?r relasjon:En differensialligning for reservoar $\mathbf{A}$
Igjen, med utgangspunkt iEn differensialligning for reservoar $\mathbf{A}$
Nei, blant annet fordi nedb?ren
En numerisk l?sning av diffligningen for reservoar $\mathbf{A}$
Vi begynner som vanlig med ? diskretisere variablene, i dette tilfellet tiden $$t \rightarrow t_i = i \cdot \Delta t$$Dette gir
En numerisk l?sning av diffligningen for reservoar $\mathbf{A}$
Til slutt l?ser vi den ukjente $A(t_{n+1})$ algebraisk og finner:Eller for ? ikke skrive for mye
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/bt.png)
Reservoar $\mathbf{B}$
Vi beskriver n? modellen for reservoar $\mathbf{B}$: Vannstand som funksjon av tid: $B(t)$ m?les i $mm$. Vann str?mmer inn i reservoaret kun gjennom bekken $Q_A(t)$. Vann str?mmer ut fra reservoaret gjennom evaporasjon $E_B$ og en bekk $Q_B$.![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/bt.png)
Reservoar $\mathbf{B}$
Regler for reservoar $\mathbf{B}$:- Innstr?mming fra $Q_A$ avtar n?r vannstanden i $\mathbf{B}$ tiltar. ($A_{sig}$ er en funksjon av $B(t)$ )
- Ingen utstr?mming til bekk n?r vannstanden er under $B_{min}$.
- $\mathbf{B}$ har i tillegg en maksgrense $B_{max}$, hvor alt overskytende vann g?r direkte i bekken $Q_{over}$
- Evapotranspirasjon er en funksjon av temperatur $E_B := E_B(T)$
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/bt.png)
Reservoar $\mathbf{B}$
Vannstanden i $\mathbf{B}$ p?virkes positivt (?kes) av innstr?mming fra $\mathbf{A}$ (gjennom bekken $Q_A$, og negativt av evaporasjon $E_B$ og utstr?mming til bekken $Q_B$ og $Q_{over}$. Alts? er $B := B(t, P_A, E_A, Q_A, Q_{over})$.![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/bt.png)
- $A_{sig}$
- $Q_{over}$
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/asig_im.png)
\begin{equation} A_{sig} = \left\{\begin{array} 1 & & B \leq B_\text{min} \\ 1- 0.25\cdot\frac{B-B_\text{min}}{B_\text{max}-B_\text{min}} & & B_\text{min} < B \leq B_\text{max} \\ 0 & & B_\text{max} < B \end{array} \right. \end{equation}
![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/bt.png)
Utstr?mming fra $\mathbf{B}$: $Q_B$
Som for $\mathbf{A}$ er ogs? str?mmen $Q_B$ ut fra $\mathbf{B}$ en funksjon av vannstanden:![](/studier/emner/matnat/ifi/IN-KJM1900/h18/forelesninger/fig/bt.png)
Utstr?mming fra $\mathbf{B}$: $Q_{over}$
I tillegg vil, dersom vannstanden overskrider et maksniv? $B_{max}$, alt overskytende vann g? rett over i bekken:Differensialligning for $B(t)$
Fra uttrykketHva skjer her?
Vi ser p? de to differensialligningene vi har funnet for modellen s? langt:Koblede ligninger
Vi ser p? de to differensialligningene vi har funnet for modellen s? langt:Numerisk l?sning av systemet
Vi kan l?se numerisk for $B(t)$ p? samme m?te som for $A(t)$, og finner at ligningssettet v?rt er\begin{align} A_{n+1} &= A_n + \cdot P_n \Delta t - \cdot E_{A,n}\Delta t - \cdot Q_{A,n} \Delta t. \\ B_{n+1} &= B_n + \cdot A_{\text{sig},n} \cdot Q_{A,n} \Delta t - \cdot E_{B,n}\Delta t - \cdot Q_{B,n} \Delta t - \cdot Q_{\text{over},n} \Delta t . \end{align}
Numerisk l?sning av ligningssettet
Hvilke st?rrelser kjenner vi ved f?rste iterasjon?\begin{align} A_{1} &= A_0 + P_0 \Delta t - E_{A,0}\Delta t - Q_{A,0} \Delta t. \\ B_{1} &= B_0 + A_{\text{sig},0} Q_{A,0} \Delta t - E_{B,0}\Delta t - Q_{B,0} \Delta t - Q_{\text{over},0} \Delta t . \end{align}
Str?m i bekken
Ved ? l?se systemet iterativt (med Eulers metode) bestemmes blant annet $A(t)$, $B(t)$, $Q_A(t)$, $Q_B(t)$ og $Q_{over}(t)$. Str?mmene inn og ut av bekken kan benyttes til ? bestemme den totale str?mmen i bekken:Kjemisk delmodell
Viktige antagelser i den kjemiske delmodellen:- Sulfat kommer med nedb?ren og f?lger vannstr?mmen uten ? reagere med omgivelsene.
- Unntak: sulfat forlater ikke systemet gjennom evapotranspirasjon.
- Elektron?ytralitet: summen av ladningene i vannet er lik 0.
- Likevekt: for andre ioner enn sulfat antar vi at de er fullstendig bestemt av likevektsreaksjoner.
Sulfattransport¶
I kolonne 1 i birkenes.dat
oppgis sulfatkonsentrasjonen i nedb?ren:
- Konsentrasjonen varierer fra dag til dag
- Av tekniske ?rsaker er den satt til negativ verdi n?r nedb?ren P = 0.
- Tenk over: kan dette bli et problem?
- Konsentrasjon $C$ angis med enheten $mol \cdot l^{-1}$ (mol per liter).
- Mengden $M$ angis (mot normalt) som $mol \cdot m^{-2}$ (mol per kvadratmeter).
Differensialligninger for sulfattransport
Ettersom evapotranspirasjon ikke transporterer sulfat ut av systemet, finner vi differensialligningene \begin{align} \frac{d}{dt} M_A & = P\cdot C_P - Q_{A}\cdot C_A, \\ \frac{d}{dt} M_B & = A_\text{sig}\cdot Q_{A}\cdot C_A - (Q_{B} + Q_\text{over})\cdot C_B, \end{align}St?rrelsene $Q_A$, $Q_B$, $Q_{over}$ og $A_{sig}$ har vi l?st for i den hydrologiske modellen, s? de resterende st?rrelsene kan l?ses numerisk.
G?r enhetene opp?
Vi ser p? enhetene i ligningen: \begin{align} \frac{d}{dt} M_A & = P\cdot C_P - Q_{A}\cdot C_A \end{align}Ja:
\begin{align} \rightarrow \frac{mol}{d?gn} = \frac{l}{d?gn} \frac{mol}{l} - \frac{l}{d?gn} \frac{mol}{l} \end{align}Den andre ligningen f?lger samme m?nster.
Differensialligninger for sulfattransport
Vi kan l?se ligningene for sulfattransport med Eulers metode, og finner da konsentrasjon og stoffmengde sulfat i reservoarene $\mathbf{A}$ og $\mathbf{B}$ og dessuten i avrenningen $Q$. Utfordringen er n? ? bestemme konsentrasjonen av den andre ionene.Elektron?ytralitet
I alle vannmengder i modellen gjelder prinsippet om elektron?ytralitet: \begin{equation} \label{elnoy} {[H^{+}] + 2[Ca^{2+}] + 3[Al^{3+}] = 2[SO^{2-}_{4}] + [HCO^{-}_{3}]} \end{equation}Her kjenner vi ved f?rste ?yekast kun sulfatkonsentrasjonen $[SO^{2-}_{4}]$.
Hvordan finnes de andre?Likevektrelasjoner:
\begin{align} K_{HCa} & = \frac{[{H^+}]^{2}}{[{Ca^{2+}}]} \\ K_{AlH} & = \frac{[{Al^{3+}}]}{[{H^+}]^{3}} \label{eq2} \\ K_{H} & = [{H^+}][{HCO^{-}_{3}}]\label{eq3} \end{align}
Elektron?ytralitet
Fra likevektene kan vi l?se for $[H^+]$ og finne \begin{align} [{H^+}] + 2\frac{[{H^+}]^2}{K_{HCa}} + 3K_{AlH}[{H^+}]^3 = 2[{SO^{2-}_{4}}] + \frac{K_{H}}{[{H^+}]}\label{elnoy2} \end{align} Vi lar $[{H^+}] := x$ og $[{SO^{2-}_{4}}] := s$ slik at uttrykket kan skrives $$3K_{AlH}x^3 + 2K_{HCa}^{-1}x^2 + x - 2s - K_{H} x^{-1} = 0$$Elektron?ytralitet
For et gitt tidssteg og reservoar kjenner vi alle st?rrelsene med unntak av $x$: \begin{align} 3K_{AlH}x^3 + 2K_{HCa}^{-1}x^2 + x - 2s - K_{H} x^{-1} = 0 \end{align} Ligningen kan derfor l?ses med Newtons metode.Ioner i vannmassene¶
Om vi kan bestemme $[{H^+}]$, s? kan vi benytte likevektene til ? finne konsentrasjonen av de andre ionene:
\begin{align} K_{HCa} & = \frac{[{H^+}]^{2}}{[{Ca^{2+}}]} \\ K_{AlH} & = \frac{[{Al^{3+}}]}{[{H^+}]^{3}} \\ K_{H} & = [{H^+}][{HCO^{-}_{3}}] \end{align}
Litt om neste uke
Vi har n? gjennomg?tt de mest sentrale delene av pensum. Det gjenst?r n? ? sette det inn i en helhetlig sammenheng, blant annet ved ? betrakte- Newtons metode og Eulers metode som rekurrensrelasjoner
- Eulers metode som numerisk integrasjon