18  Der Wasserstand eines Sees

Angenommen, du sitzt nach deinem Abschluss in einem Umweltplanungsbüro und erhältst an deinem ersten Tag einen Anruf eines Gemeindebürgermeisters aus den Voralpen. Er berichtet, dass es aufgrund des Klimawandels immer mal wieder zur Wasserknappheit im örtlichen See kommt. Er erzählt außerdem, dass zukünftig eine Industrieanlage geplant ist und er Sorge hat, dass diese Schadstoffe in den See einleiten könnte. Der See ist für die Trinkwasserversorgung von großer Bedeutung und er möchte wissen, ob du helfen kannst. Nun erwartet er vermutlich kein Modell, sondern Vorhersagen und bestenfalls Handlungsanweisungen. Da du nicht weißt, was seitens eures Büros so üblich ist und du dem Bürgermeister einen Rückruf eures Chefs angekündigt habt, macht du dir also den Rest des Tages Gedanken, wie du helfen kannst. Viel hast du bisher noch nicht erfahren, aber du findest per KI ein paar Eckdaten zu diesem See heraus. Viel wichtiger ist aber zunächst eine Zielfrage. Vor der Vorhersage und darauf aufbauenden Handlungsanweisungen, muss zunächst der Status Quo stehen. Was du brauchst ist der Wasserstand des Sees im Jahresverlauf. Hier reicht es aber nicht, dass du dazu Messdaten vom örtlichen Wasserwerk bekommst. Dur musst die Ursache für die Wasserstandsschwankungen richtig beschreiben.

Abbildung 18.1: Der Modellsee

Die Wasserstandsänderungen sollen von Tag zu Tag berechnet werden. Die Oberfläche des (kreisförmigen) Sees betrage bei mittlerem Wasserstand \(2\,km^2\). Seine normale Tiefe beträgt \(6\,m\). Der See hat einen Zufluss und einen Abfluss. Der Zufluss ist in der kälteren, regenreichen Zeit hoch, insbesondere gibt es ein Frühjahrshochwasser (nach der Schneeschmelze). Im Sommer kann es eine Periode ohne Zufluss geben. Der Abfluss hängt vom Wasserstand ab. Sinkt dieser unter eine Schwelle, versiegt er ganz. Durch die temperaturabhängige Oberflächenverdunstung treten im Sommer Wasserverluste auf.

18.1 Verbale Beschreibung und Ziel des Modells

Eine verbale Beschreibung des Modells ist dann unerlässlich, wenn man nicht allein am Modell arbeitet, oder wenn der Modellierer einen Auftraggeber hat. Gegebenenfalls müssen Schnittstellen zu anderen Modellen definiert werden. Die verbale Beschreibung ist eine Vorstufe des mathematischen Modells: ein konzeptionelles Modell. Insbesondere ist es wichtig, die Ziele und die Fragen, die mit dem Modell beantwortet werden sollen, gut zu formulieren. Man modelliert ja nicht einfach “ein System”, sondern modelliert es für einen Zweck, der die Art und den Gültigkeitsbereich des Modells bestimmt.

18.2 Festlegung der Zustandsvariablen und Schematisierung des Modellsystems

Dynamische Systeme haben einen zeitabhängigen “Zustand”. Der Zustandsraum ist die Menge aller möglichen Zustände, die das Modellsystem im Prinzip annehmen kann. Der aktuelle Zustand des Systems zu einem ausgewählten Zeitpunkt t ist ein Punkt im Zustandsraum. Die Zustandsvariablen (auch Systemvariablen oder prognostische Variablen) beschreiben also die zeitliche Entwicklung des Systems. Man kann sie meistens nicht direkt angeben, sondern stellt Gleichungen auf, die die Änderung des Zustands angeben. Hierzu werden die Prozesse bestimmt, die die Zustandsvariablen ändern.

Im Falle des Sees bietet sich auf den ersten Blick als einzige Zustandsvariable der Wasserstand an. Später kann als zweite Zustandsvariable die Schadstoffkonzentration hinzukommen. Als Pegel oder Wasserstand \(W\) bezeichnen wir die aktuellen Tiefe des Sees¹. Es könnte aber auch ein Pegel sein, der sich durch eine Nullpunktverschiebung von der aktuellen Tiefe des Sees unterscheidet. Unser Ziel ist es, den jahreszeitlichen Verlauf des Wasserstands zu ermitteln - wir suchen W(t).

18.3 Festlegung des Modellsystems - Topographie des Sees

Wir nehmen vereinfachend an, dass der See einen parabelförmigen vertikalen Querschnitt hat, der See sei also ein Paraboloid. Pegelnull wird hier vereinfachend als 0 angenommen, der Pegel steht sozusagen im Scheitelpunkt des Sees. Im allgemeinen ist “absolut” Null nicht dasselbe wie Pegelnull.

Abbildung 18.2: Schematische Darstellung des Modellsees mit parabelförmigen Querschnitt. \(W_N\) gibt den Normwasserstand an. Bei diesem Wasserstand hat der See den Radius \(R_N\) also die Fläche \(F_N=\pi R_N^2\). Der rote Punkt markiert die Höhe des Abflusses.

Mit dem Wasserstand verknüpft ist die Seeoberfläche F, die zur Bestimmung der Verdunstung benötigt wird. Wir brauchen also die Seeoberfläche \(F\) als Funktion des Wasserstandes \(W\). Zuerst betrachten wir die Fläche als Funktion des Radius:

\[F(R)=\pi R^2\]

Wir nehmen an, dass ein Referenz-Wasserstand \(W_N\) und \(F_N\) die zugehörige Fläche (z.B. aus Satellitendaten) bekannt sind

\[ W(F_N)=W_N \text{ und } F_N=\pi R_N^2 \tag{18.1}\]

Man kann nun den Wasserstand \(W\) als Funktion des Radius \(R\) mit dem Ansatz

\[W(R) = a \cdot R^2\]

beschreiben. Die Konstante \(a\) ergibt sich aus der Referenzbeziehung -@ eq-see:ref:

\[W(R_N)=a \cdot R_N^2 =W_N\quad \Rightarrow \quad a=\frac{W_N}{R_N^2} = \pi \frac{W_N}{F_N}\, ,\]

so dass

\[W(R) = \frac{W_N}{R_N^2} \cdot R^2 = \frac{W_N}{F_N}\cdot F =W(F)\]

Man erhält die Fläche in Abhängigkeit des Wasserstands durch Auflösen nach F:

\[ F\left(W\right)=F_N \cdot \frac{W}{W_N} \tag{18.2}\]

Die Fläche F ist also proportional zur Tiefe \({W}\).

18.4 Festlegung der Prozesse

Die Prozesse sind die Vorgänge, die die Zustandsvariablen ändern. Der Wasserstand \(W\) wird durch die Prozesse Zufluss, Abfluss und Verdunstung geändert:

Abbildung 18.3: Schematische Darstellung des Modellsees. Die Pfeile stellen hier materielle Flüsse dar (Wassermengen pro Zeit).

Alle drei Größen haben die Dimension Volumen pro Zeiteinheit, also z.B. die Einheit \(m^3d^{-1}\) (Kubikmeter pro Tag). Es ist wichtig, sich immer über Dimension und Einheiten aller auftretender Größen im Klaren zu sein und bei jedem Term, bei jeder Gleichung sofort die Einheiten auf Konsistenz zu überprüfen! Das ist ein wichtiges Hilfsmittel, um Fehler beim Aufstellen der Gleichungen frühzeitig aufzuspüren.

Die Modellgleichungen verbinden die Zustandsvariablen mit den Prozessen in Form von “dynamischen” Gleichungen, das sind Gleichungen, die die Änderungen der Zustandsvariablen beschreiben. In der Regel sind das Differentialgleichungen oder Differenzengleichungen. Die Prozesse ändern den Wasserstand, indem sie das Volumen ändern:

In unserem Modell soll die Zeit in Tagesschritten \(\Delta t\) fortschreiten (Differenzengleichungsmodell):

\[ \Delta V=\left(zu - ab -vd \right)\cdot\Delta t \tag{18.3}\]

Die Multiplikation mit \(\Delta t\) auf der rechten Seite ist wichtig: Zum einen werden die Prozesse mit Dimension Volumen pro Zeit durch die Multiplikation mit einer Zeit in Volumen verwandelt, wie es die linke Seite erfordert. Zum anderen erkennt man eine “Linearität”: Die Volumenänderung ist proportional zur Länge des Zeitschrittes. Genau hierin besteht die Approximation der Differenzengleichung, denn offensichtlich ist diese Linearität nur für kleine Zeitschritte \(\Delta t\) akzeptabel, nämlich solange sich \(zu\), \(ab\) und \(vd\) in \(\Delta t\) nur unwesentlich ändern. In Abbildung 18.3 ist die dynamische Gleichung in Form eines Diagramms dargestellt.

Sprachlich unscharf und sogar falsch wäre die Ausdrucksweise: “Das Volumen hängt von Zufluss, Abfluss und Verdunstung ab”. Richtig ist aber: “Das Volumen wird durch Zufluss, Abfluss und Verdunstung geändert”. Das Volumen \(V\) ist keinesfalls eine “Funktion” von \(zu\), \(ab\) und \(vd\). Dann müsste sich nämlich aus den augenblicklichen Werten von \(zu\), \(ab\) und \(vd\) auf das augenblickliche Volumen schließen lassen, was nicht der Fall ist. Man sieht es schon daran, dass das Volumen unterschiedlich sein kann, wenn z.B. alle drei Flüsse verschwinden. Nicht \(V\) selbst sondern die zeitliche Änderung des Volumens ist eine “Funktion” der Prozesse!
Als Zustandsvariable hatten wir aber den Wasserstand \(W\) und nicht das Volumen gewählt!

Abbildung 18.4: Im Zeitintervall \(\Delta t\), z.B. in einem Tag, ändert sich das Volumen um
\(\Delta V =\left(zu -ab -vd \right)\cdot \Delta t\). Bei gegebener Oberfläche \(F\) des Sees sind Volumenänderung und Wasserstandsänderung durch \(\Delta V=F\cdot \Delta W\) verknüpft. Damit ergibt sich Gleichung 18.4.

Der Zusammenhang zwischen \(V\) und \(W\) über die Oberfläche \(F\) ist annähernd (diskret) durch

\[\Delta V = F\cdot \Delta W\]

gegeben. Damit lässt sich Gleichung 18.3 entsprechend umschreiben:

\[ \Delta $W$ =\frac{1}{F(W)}(zu-ab-vd)\cdot \Delta t \tag{18.4}\]

Dies ist die gesuchte Modellgleichung. Die Prozesse \(zu\), \(ab\) und \(vd\) müssen jetzt noch präzisiert werden.

Bemerkung: Der integrale Zusammenhang zwischen \(V\) und \(W\) über die Oberfläche \(F\) ist durch \(V\left(W\right) =\int\limits^W_0 F\left(w\right)\, dw\) gegeben. Er hat die differentielle Form \(\dot V =\frac{dV}{dW}\cdot \frac{dW}{dt}= F\cdot \dot W\).

18.5 Festlegung der Randbedingungen und Antriebe (Forcing)

Bei einem “System” muss man immer klar trennen, was dazu gehört und was als extern anzusehen ist. Die Umgebung beeinflusst die Entwicklung des Modellsystems, aber – und das ist eine Folge der (künstlichen) Modellabgrenzung – das Modellsystem kann das “Externe” nicht ändern, auch wenn es “in Wirklichkeit” eine Rückwirkung über die Systemgrenzen hinaus gibt. Soll diese modelliert werden, sind die Systemgrenzen weiter zu fassen. Die externen Bedingungen, die für den Fortgang der internen Entwicklung wesentlich sind, werden “Antriebskräfte” genannt (“driving forces”). Bei Ökosystemen sind wichtige Antriebskräfte Licht, Temperatur und Niederschlag. Sie ändern sich mit der Jahreszeit und der Tageszeit. Aber auch andere “Randbedingungen” gehören zu den Antriebskräften. Das sind z.B. Konzentrationen von Stoffen außerhalb der räumlichen Systemgrenzen, wenn sie diese durch Transportprozesse überschreiten und dadurch im System die Zustandsvariablen ändern.

Als Forcing werden also alle Größen bezeichnet, die das System von “außen” beeinflussen. Im Falle des Sees sind das die Niederschläge, die den Zufluss bestimmen, und die Temperatur, die die Verdunstung reguliert. Beide ändern sich mit der Jahreszeit. Der Abfluss hängt im Gegensatz zum Zufluss vom Systemzustand (Wasserstand) selbst ab. Die Verdunstung ist neben der Temperatur-Abhängigkeit proportional zur Oberfläche. In der Realität wirkt ein großer See durch Verdunstung auf das Regionalklima und somit auf Wolken, Niederschlag, Zufluss, Luftfeuchtigkeit und Verdunstung zurück. Durch Festlegung der Systemgrenzen wird diese Rückwirkung im Modell aber gerade ausgeschlossen.

18.5.1 Niederschlag und Zufluss

Der Niederschlag bestimmt letztlich den Zufluss. Will man aus Niederschlagswerten, die meistens in \(mm\,d^{-1}\) angegeben werden, den Zufluss bestimmen, muss man das Einzugsgebiet des Sees kennen und die Zeitverzögerung, mit der der Niederschlag tatsächlich im See ankommt. Statt des Niederschlags geben wir daher vereinfachend den Zufluss selbst als Antrieb vor. Der Zufluss ändert sich natürlich von Tag zu Tag, häufig bekommt man als Daten aber nur monatliche Mittelwerte. Diese müssen dann in Tageswerte umgerechnet werden. Um ein Modell zu testen, nimmt man auch sogenannte klimatologische Daten (Mittel über 30 Jahre). Diese Daten werden dann bei der Simulation von Jahr zu Jahr wiederholt. Später kann man diese Daten dann durch die aktuellen Werte ersetzen.

Abbildung 18.5: Monatsmittelwerte des Zuflusses.

Die Dimension der Zuflussdaten \(v_{zu}\) ist Volumen pro Zeit. In Abb. fig-see:2 sind die Daten in der Einheit \(m^3\,d^{-1}\) angegeben.

18.5.2 Temperatur

Für die Temperatur geben wir die mittlere, jährliche Temperatur \(T_m\), die Amplitude der Temperaturschwankung \(T_a\) und den kältesten Tag \(t_k\) im Jahr vor und bestimmen die Temperaturkurve wie in Abschnitt 2.1. Später kann man auch diese Daten durch reale Daten ersetzen.

\[ T\left(t\right)=T_m-T_a\cdot \cos \left(2\pi \frac{t-t_k}{360\,d}\right) \tag{18.5}\]

Abbildung 18.6: Temperaturverlauf im Jahresgang in ° C. Kältester Tag: 1. März.

18.6 Prozessbeschreibungen

18.6.1 Zufluss

Der Zufluss ist direkt durch die Daten gegeben (Abbildung 18.5). Modelltechnisch muss man sicherstellen, dass der Zufluss am Ende eines Monats auf den Wert des Folgemonats springt. Man muss also zu einem gegebenen Tag \(t\) den Monat \(m\) ausrechnen. Da wir von einer Monatslänge von 30 Tagen ausgehen, ist dies recht einfach. Bei realistischer Betrachtung müsste man die unterschiedlichen Monatslängen sowie Schaltjahre berücksichtigen. Hier kann man es einfach durch die Modulo-Funktion lösen²:

\[m=\left [\frac{1}{30}\left((t-1) \mod 360 \right)+1\right ] \nonumber\]

Der Zufluss beträgt somit einfach

\[ zu =v_{zu}\left(m\right) \tag{18.6}\]

wobei \(v_{zu}\) der mittlere monatliche Zufluss ist. Der Zufluss hängt nicht von der Zustandsvariablen \(W\) ab.

18.6.2 Abfluss

Wenn der Wasserstand unter die Abflusshöhe sinkt, versiegt der Abfluss. Je höher der Wasserstand ist, desto schneller fließt das Wasser ab. Nach dem Gesetz von Toricelli ist die Abflussgeschwindigkeit proportional zur Wurzel aus der Höhe der Wassersäule über dem Abfluss. Man kann folgenden einfachen Ansatz wählen:

\[ ab\left(W\right) = F_{ab} \cdot \sqrt{max(0, 2 g \cdot(W-W_A))} \cdot sec_{day} \tag{18.7}\]

wobei \(g \approx 9.81 m\,s^{-2}\) die Gravitationskonstante, \(F_{ab}\) (in \(m^2\)) eine Proportionalitätskonstante, die von der Querschnittsfläche des Abflusses abhängt, ist und \(sec_{day}=86400\) (in \(s\,d^{-1}\)) die Anzahl Sekunden in einem Tag angibt. Letztere Größe stellt sicher, dass der Abfluss die Einheit \(m^3\,d^{-1}\) hat.

In der Realität spielen viele weitere Faktoren, wie die Beschaffenheit des Abflusses, die Viskosität des Wasser, Turbulenzen etc., eine Rolle. Die Messung von Abfluss- bzw. Durchflussraten ist messtechnisch nicht einfach, daher wird man nur selten hochaufgelöste Daten für einen gewünschten Zeitraum finden.

18.6.3 Verdunstung

Die Verdunstung \(vd\) ist proportional zur aktuellen Oberfläche \(F\) des Sees und hängt somit vom aktuellem Wasserstand \(W\) ab. Sie steigt mit der Temperatur an und hängt in der Realität außerdem vom Wind und von der Luftfeuchtigkeit ab. Wir ignorieren diese zusätzlichen Abhängigkeiten.

Wir nehmen an, dass die Abnahme des Wasserstands durch Verdunstung bei \(20 ° C\) bekannt ist. Diese Verdunstungskonstante \(vd_{20}\) hat die Dimension Länge pro Zeit. Verdunstung wird meistens in mm angeben, so dass \(vd_{20}\) die Einheit \(mm\, d^{-1}\) hat. Um nun von der Verdunstungskonstante auf das Verdunstungsvolumen zu kommen, muss diese mit der aktuellen Seeoberfläche multipliziert werden. Damit ist \(vdmm_{20}\cdot F\) das Volumen, das bei \(20 ° C\) pro Tag verdunstet. Diese Größe hat bisher allerdings die Einheit \(mm\, m^2 \, d^{-1}\) und muss daher noch durch 1000 geteilt werden.

Um nun die Temperaturabhängigkeit zu berücksichtigen, benötigen wir eine (dimensionslose) Funktion der Temperatur, die gerade bei \(20 ° C\) den Wert 1 annimmt und die Verdunstung physikalisch beschreibt. Durch Multiplikation mit dieser Funktion, nennen wir sie \(e_T\), erhalten wir die gesuchte Prozessbeschreibung für die Verdunstung

\[ vd\left(t,W\right)=vd_{20} \cdot F\left(W\right) \cdot e_T\cdot 10^{-3} \tag{18.8}\]

Man beachte, dass \(e_T\) von der Temperatur abhängt und diese wiederum von der Zeit. Es fehlt nun noch die Präzisierung von \(e_T\). Dies ist nicht einfach, da die Verdunstung von vielen Faktoren, wie z.B. der Temperatur der Luft, der Luftfeuchtigkeit, Wind- und Wellengang abhängt. Es fehlt nun noch die Präzisierung von \(e_T\). Dies ist nicht einfach, da die Verdunstung von vielen Faktoren, wie zB. dem Dampfdruck abhängt. Vereinfachend gehen wir von einem exponentiellen Zusammenhang aus, indem wir annehmen, dass sich die Verdunstung bei einer Temperaturerhöhung um \(T_{vd}\) gerade verdoppelt:

\[ e_T= 2^\frac{T-20 ° C}{T_{vd}} \tag{18.9}\]

So ist \(e_T\) dimensionslos und nimmt bei \(20 ° C\) den Wert 1 an. Hinter diesem Ansatz steckt ein wichtiges Prinzip, das man beim Aufstellen mathematischer Modelle beachten sollte: Die Formeln werden so strukturiert, dass dimensionsbehaftete Größen durch dimensionslose Faktoren modifiziert werden. So kann man die Modifikation ggf. anpassen, wenn man neue Erkenntnisse gewonnen hat, ohne dass man sich ein weiteres Mal Gedanken um die Einheiten machen muss.

18.7 Bestimmung der Parameterwerte

Den Modellparametern, müssen Werte zugeordnet werden. Viele Werte kann man aus der Literatur, oder annähernd aus Lehrbüchern und durch Expertenbefragung gewinnen. Einige Parameterwerte muss man schätzen, da sie nicht bekannt sind. Diesen “Schätzen” von Parameterwerten, nennt man “educated guessing”. Es handelt sich nicht um ein blindes Herumraten, sondern um die Umsetzung des Systemverständnisses, das sich auf den bisherigen Schritten und durch Literaturstudium eingestellt hat. Regelmäßig ist ein Teil der Parameter in den Modellgleichungen erst durch Aufstellung der Gleichungen definiert worden. Für diese Parameter können vorab keine Werte durch Messungen bestimmt worden sein! Durch die Modellierung liegt nun eine theoretische Vorstellung des Systems vor, die es erlaubt, Experimente zu planen, mit denen diese neuen Parameter gemessen werden könnten. Zunächst aber muss man diese Werte schätzen. In Tabelle 18.1 sind alle Größen, die den See beschreiben, Parameter und Antriebe mit Werten, Einheiten und Bedeutung zusammengefasst.

Symbole	Einheiten	Werte	Erklärung
Konstanten
\(g\)	\(mm \cdot m\, s^{-2}\)	9.81	Gravitationskonstante
\(sec_{day}\)	\(s\,d^{-1}\)	86400	Sekunden eines Tages
Seegeometrie
\(W_N\)	\(m\)	6	Normpegel
\(F_N\)	\(m^2\)	2\(\cdot 10^6\)	Seeoberfläche bei Normpegel
\(W_A\)	\(m\)	2	Lage des Abfluss
Dynamische Parameter
\(F_{ab}\)	\(m^2\)	0.04	“Abflussquerschnitt”
\(vd_{20}\)	\(mm \cdot d^{-1}\)	15.0	Verdunstung bei 20° C
\(T_{vd}\)	\(° C\)	5.0	Verdopplung der Verdunstung
Antriebskräfte
\(v_{zu}\)	\(10^3\,m^3\cdot d^{-1}\)	40 20 140 60 20 10 0 0 0 20 50 40	Monatsmittel des Zuflusses (Jan - Dez)
\(T_m\)	\(° C\)	15	Mittlere Jahrestemperatur
\(T_a\)	\(° C\)	12	Temperaturamplitude
\(t_k\)	\(d\)	60	Kältester Tag (1.3.)

Tabelle 18.1: Zusammenfassung der Seegrößen, Parameter und Antriebe mit Einheiten, Werten und Bedeutung

18.8 Prüfung des Modells auf Konsistenz

Der See hat bei dem Wasserstand \(W_N=6\,m\) eine Oberfläche von \(F_N=2\,km^2\) . Bei kreisrunder Seefläche und parabelförmigen Tiefenprofil ergibt das einNormseevolumen³ von \(6\cdot 10^6\,m^3\). Der Zufluss ist durch die Monatsdaten vorgegeben. Der jährliche Zufluss beträgt danach \(30 \cdot \sum \limits_{m=1}^{12} v_{zu}(m) = 12 \cdot 10^6 \; m^3\). Dies bedeutet, dass der Zufluss das Zweifache des Normseevolumens beträgt.

Bei Normpegel beträgt der Abfluss nach Gleichung 18.7 in etwa \(3\cdot 10^4\,m^3\, d^{-1}\). Dies ergibt ein geschätztes jährliches Abflussvolumen von \(11 \cdot 10^6\,m^3\) was in etwa dem Zufluss entspricht.

Die Verdunstung konzentriert sich auf die vier Sommermonate, in denen die Temperaturen um \(25° C\) erreicht werden. In diesen vier Monaten (\(120\,d\)) ergibt sich nach Gleichung 18.8 bei halber Normfläche eine jährliche Verdunstung von in etwa \(3.6 \cdot 10^6\,m^3\). Die Verdunstung ist also geringer als Zu- und Abfluss.

Zufluss, Abfluss und Verdunstung werden sich über das Jahr betrachtet nicht ganz ausgleichen. Die Größenordnungen aller drei Flüsse erscheinen aber einigermaßen plausibel. Damit sind die gewählten Parameter eine vernünftige Basis für erste Simulationen.

18.9 Programmierung, Simulation, Sensitivitätsuntersuchung

Nachdem die drei Prozesse durch Gleichungen präzisiert wurden, lässt sich das Modell als Ganzes durch Zusammenfassung der Gleichungen darstellen. Dies ist die mathematische Formulierung des Modells, das “Entwicklungsgesetz”, da es den Zustand zur Zeit \(t+\Delta t\) aus dem Zustand zur Zeit \(t\) zu berechnen erlaubt. Aus gegebenen "‘alten"’ Werten werden “neue” berechnet, und das wird iterativ wiederholt. Schematisch sieht das wie folgt aus (Pseudocode):

W(1) =5.0 Wasserstand zur Zeit 0
dt =1.0   Modellzeitschritt    
     
Schleife ueber die Zeitschritte i
      t = i*dt 
      F =...
      T =...
      ...
      zu = ...
      ab =...
      vd =...
      Delta_W = (zu-ab-vd)/F 
      W(i+1) = W(i)+Delta_W*dt   
Ende der Schleife

Bei einer Iteration muss also ein “Iterationsanfang” vorgegeben werden. Zu Beginn muss für einen Anfangszeitpunkt \(t_0\) der zugehörige Wasserstand \(W_0\) vorgegeben werden. Dann beginnt die Iteration (Schleife): Für gegebene Werte t und \(W(t)\) werden zuerst die Hilfsgrößen Fläche \(F\) und Temperatur \(T\) und Monat \(m\) berechnet. Diese Hilfsgrößen werden für die Berechnung der Prozesse zu, ab und vd benötigt. Zum Schluss wird mit ihnen die Änderung des Wasserstandes \(\Delta W\) und der neue Wasserstand \(W\left(t+\Delta t\right)=W+\Delta W\) bestimmt. Die Zeit schreitet um \(\Delta t\) fort. Nach dem Durchlaufen kann diese Abfolge mit den “neuen” Werten t und \(W\) erneut gestartet und immer wiederholt werden, bis ein Endzeitpunkt erreicht ist. Der Iterationsanfang bei \(t_0\) mit dem Startwert \(W(t_0)=W_0\) legt so die Iteration fest.

Man kann davon ausgehen, dass die Modellresultate der ersten Versuche nicht das gewünschte Ergebnis zeigen. Nun kommt es darauf an, herauszufinden, woran das liegt. Ist es ein Programmierfehler, ein Fehler in den Parameterabschätzungen oder ein Denkfehler beim Aufstellen der Modellgleichungen? Wenn es sich um Programmierfehler handelt, so erfordert es Geduld, Übersicht und detektivische Kleinarbeit, sie in einem längeren Code aufzuspüren. Fehlersuche ist anspruchsvoll und wird durch Erfahrung erleichtert.

In Abbildung 18.7 sind die Simulationsergebnisse für verschiedene Startwerte \((t_0,H_0)\) (Zeit und Wasserstand) dargestellt.

Abbildung 18.7: Wasserstand des Sees im Jahresgang. Unterschiedliche Anfangspegel “vergisst” das System nach wenigen Monaten. Es stellt sich bei periodischen Antrieben (klimatologisches Jahr) nach einiger Zeit eine periodische Lösung ein.

Man gewinnt aus Abbildung 18.7 eine wesentliche Erkenntnis: Das System “vergisst” den Anfangswert nach einiger Zeit, hier bald nach dem ersten Hochwasser. Solche Systeme nennt man “dissipativ”. Nach einiger Zeit stellt sich in jedem Fall die gleiche periodische Lösung ein, da die Antriebe (Temperatur und Zufluss) von Jahr zu Jahr wiederholt werden (klimatologische Mittelwerte). Es gibt offenbar eine ausgezeichnete periodische Lösung, die zum Anfangswert \(W_0\approx2.9\,m\) gehört, und zu der alle anderen Lösungen hinstreben.

In Abbildung 18.8 sind die Simulationsergebnisse für das dritte Simulationsjahr und die zugehörigen Flüsse dargestellt. Man erkennt sofort den Einfluss des Zufluss auf den Wasserstand. Im August und September sind Zufluss und Abfluss auf Null. Die Verdunstung spielt kaum eine Rolle.

Abbildung 18.8: Oben: Wasserstand des Sees des 3. Jahres. Die schwarze Linie markiert die Höhe der Abflussschwelle. Unten: Zufluss, Abfluss und Verdunstung des 3. Jahres.

Wenn das Modell befriedigend läuft, kann man durch Parametervariation untersuchen, wie sensitiv das Modell reagiert. Durch Sensitivitätsanalysen sind die Parameter zu identifizieren, die das Simulationsergebnis empfindlich beeinflussen. Dazu verändert man einzelne Parameter geringfügig nach oben und unten und vergleicht die zugehörigen Simulationen miteinander. Gerade bei geschätzten Parametern ist dies wichtig, um ein Gefühl dafür zu bekommen, wie “wichtig” ein unsicherer Parameter ist. Insbesondere möchte man wissen, in welchem Bereich sich die Parameterwerte bewegen, die zu sinnvollen Ergebnissen führen. Ein Verständnis in die Plausibilität des Modells erhält man auch dadurch, dass man die Grenzen des Modells auslotet.

18.10 Modellverbesserung

Spätestens an dieser Stelle erweist sich der Modellierungsprozess als Iterationsschleife. Man muss wieder und wieder die Modellgleichungen in Frage stellen (nicht nur den Code oder die Parameterwerte) und die Implementierung überprüfen. Eventuell müssen sogar die Ziele des Modells und dessen Gültigkeitsbereich hinterfragt werden. Der Vergleich zwischen Simulation und Messdaten wird zur Verbesserung des Modells herangezogen.

18.11 Anwendungen des Modells, Dokumentation der Ergebnisse

Wenn das Modell vertrauenswürdig “läuft”, kann man es in verschiedenen Anwendungen nutzen, um weitere Einsichten in das System zu erhalten. Man kann Szenarien aller Art berechnen, und insbesondere das Modell unter Bedingungen laufen lassen, die man dem realen System nicht zumuten darf. Die Dokumentation des Modells und die Publikation der Ergebnisse sind unabdingbar.

18.12 Schadstoffeinleitung

Es wird angenommen, dass ein Schadstoff in den See eingeleitet wird. Schadstoffe werden im allgemeinen in Konzentrationen (Masse pro Volumen) angegeben. Daher ist es nun sinnvoll, statt des Wasserstandes das Volumen des Sees als Zustandsvariable zu betrachten. Der Wasserstand kann dann nach Gleichung 18.4 ableitet werden:

\[\Delta W =\frac{1}{F(W)}(zu-ab-vd)\cdot \Delta t\]

geht unter Berücksichtigung der Beziehung \(\Delta V=F\cdot \Delta W\) über in

\[\Delta V =(zu-ab-vd)\cdot \Delta t\]

Weiterhin werden wir das Modell als DGL-Modell betrachten, um im folgenden Instabilitäten durch zu schnelle Änderungen zu vermeiden:

\[ \dot V=zu-ab-vd \tag{18.10}\]

Das Modell wird nun um eine Zustandsvariable erweitert. Neben dem Wasservolumen \(V\) wird die Konzentration \(C\) oder die Menge \(M\) eines Schadstoffs im See berücksichtigt. Man kann \(C\) oder \(M\) als zweite Zustandsvariable wählen. Da im allgemeinen die Menge des eingeleiteten Schadstoffs eher messbar ist und das Volumen über die Zeit veränderlich ist, ist es sinnvoll die Masse als Zustandsvariable zu wählen und die aktuelle Konzentration aus aktueller Masse und aktuellem Volumen abzuleiten (\(C=\frac{M}{V}\)). Für die Schadstoffmenge \(M\) muss also eine Modellgleichung aufgestellt werden. Die Schadstoffmenge ändert sich durch die Prozesse Einleitung \(ein\) und Ausstrom \(aus\). Die Einleitung soll vom Seeufer erfolgen und unabhängig vom Wasserzufluss in den See sein (z.B. ein Industrieabwasser als eingeleitete Menge pro Tag). Die Einleitung sei also eine Funktion der Zeit. Der Schadstoffausstrom hängt vom Wasserabfluss und der aktuellen Schadstoffkonzentration \(C\) sec-ab-

\[aus=ab\cdot C\nonumber\]

Wir nehmen weiter an, dass der Zufluss in den See schadstofffrei ist und auch keine Schadstoffverdunstung oder Ausfällung stattfindet. Die Modellgleichungen sind dann durch Gleichung 18.11 gegeben.

\[\begin{aligned} \dot V&=zu-ab- vd\\\nonumber \dot M&=ein - aus \nonumber \end{aligned} \tag{18.11}\]

Wenn die Konzentration \(C\) statt der Menge \(M\) verwendet werden soll, so ergibt sich wegen \(M=V\cdot C\) und \(\dot M=\dot V\cdot C + V\cdot\dot C\) das zu 18.11 äquivalente System 18.12.

\[\begin{aligned} \dot V&=zu-ab-vd\\\nonumber \dot C &=\frac{1}{V} \left(ein-C\cdot zu+C\cdot vd\right)\\ \nonumber \end{aligned} \tag{18.12}\]

Die Gleichungen zeigen, dass die Menge des Schadstoffs durch den Abfluss abnimmt, dass aber der Abfluss die Konzentration nicht beeinflusst. Umgekehrt wird die Menge durch Zufluss und Verdunstung nicht geändert, während die Konzentration durch Zufluss abnimmt (Verdünnung) und durch Verdunstung zunimmt (Eindickung). Wenn bei Trockenfall des Sees das Volumen V verschwindet, steigt die Konzentration C über alle Grenzen (nicht aber die Menge M). Unter der Annahme, dass zu Simulationsbeginn der See nicht verschmutzt ist und pro Tag 10 g des Schadstoffs eingeleitet werden, ergibt sich der in Abbildung 18.9 dargestellte Verlauf.

Abbildung 18.9: Schadstoffeinleitung in einen See. Simulation über 5 Jahre. Dargestellt ist der Pegelstand (oben) und die Schadstoffkonzentration (unten) für das 5. Jahr bei einer Einleitung von 10 g pro Tag. Zu Beginn der Simulation war der See schadstofffrei.

18.12.1 Vermeidung hoher Schadstoffkonzentrationen

Die Einleitung des Schadstoffs kann als Antrieb betrachtet werden, solange sie nicht direkt vom Wasserstand W oder von der Konzentration C im See abhängt. Die Einleitung kann optimiert werden, wenn wir unser Wissen über das System nutzen. Es könnte sinnvoll sein, die Einleitung auf die niederschlagsreichen Monate November bis April zu beschränken und dafür zu verdoppeln oder auf einen kurzen Zeitraum vor dem Hochwasser zu konzentrieren. In der Simulation werden drei Szenarien verglichen.

1. konstante Einleitung von 10 g pro Tag

2. Einleitung von 20 g pro Tag von November bis April

3. Einleitung von 30 g pro Tag von Dezember bis März

Die Gesamtmenge der Einleitung bleibt so erhalten.

Das Ergebnis der Simulationen ist in Abbildung 18.10 dargestellt, wobei für alle drei Szenarien der See anfangs nicht verunreinigt ist. Es zeigt sich, dass man die Konzentrationsspitzen im Sommer deutlich vermindern kann, wenn die Einleitung auf die Wintermonate konzentriert wird. Der starke Konzentrationsanstieg erfolgt im Sommer immer, auch wenn in dieser Zeit gar keine Einleitung stattfindet. Das ist die Folge der Verdunstung des Wassers, bei der der Schadstoff zurückbleibt. Im Winter nimmt die Konzentration ab, auch wenn eingeleitet wird, weil sich das Volumen durch Wasserzufluss schnell vergrößert.

Abbildung 18.10: Simulation der Schadstoffeinleitung in den See über 5 Jahre. Die gleiche Menge eines Schadstoffs wird auf 12 Monate verteilt (blau), oder auf 6 (gelb) bzw. 4 (rot) Wintermonate konzentriert. Die mittlere Einleitung beträgt jeweils 10 g pro Tag. Dargestellt ist jeweils das 5. Jahr. Zu Beginn der Simulation war der See schadstofffrei.

18.12.2 Bakterieller Abbau

Findet im See ein bakterieller Abbau des Schadstoffs statt, so verringert dieser die Schadstoffmenge:

\[\begin{aligned} \dot V&=zu-ab- vd\\\nonumber \dot M&=ein - aus -abbau\nonumber \end{aligned} \tag{18.13}\]

Wir nehmen an, dass der Abbau temperaturabhängig und proportional zur aktuellen Schadstoffmenge ist. Ist \(r\) die bakterielle Abbaurate (in \(d^{-1}\)) bei \(10 ° C\) und \(e_T\) der Temperaturfaktor nach Gleichung 18.9 so kann man denn Abbau durch

\[\begin{aligned} abbau=r\cdot e_T \cdot M\nonumber \end{aligned}\]

beschreiben. Damit ist die Temperaturabhängigkeit des bakteriellen Abbaus entsprechend der van’t Hoffschen Regel (Sommer 2005) berücksichtigt.

Abbildung 18.11: Simulation der Schadstoffeinleitung in den See über 5 Jahre bei einer kontinuierlichen Schadstoffeinleitung von 10g pro Tag. Dargestellt ist die Entwicklung der Konzentration im See ohne (blau) und mit (rot) bakteriellem Abbau von \(1 \%\) der aktuellen Schadstoffmenge im See pro Tag. Zu Beginn der Simulation war der See schadstofffrei.

Die Ergebnisse der Simulation sind in Abbildung 18.11 dargestellt. Man sieht deutlich die Wirkung des bakteriellen Abbaus. Der Spitzenwert von über \(100\,mg\, m^{-3}\) wird durch den bakteriellen Abbau von 1% pro Tag auf unter \(1\,mg\, m^{-3}\) reduziert.

Wenn es sich bei dem betrachteten Schadstoff um einen Nährstoff handelt, wird die Bakterienzahl mit zunehmender Nährstoffverfügbarkeit zunehmen. Dies wird hier nur implizit berücksichtigt, da der bakterielle Abbau als proportional zur Schadstoffmenge angenommen wird. Genauer wäre es, die Bakterien als weitere dynamische Größe (Zustandsvariable) zu betrachten und ihr Wachstum von weiteren Faktoren wie der Verfügbarkeit organischen Materials abhängig machen.

Sommer, U. 2005. Biologische Meereskunde. 2. Aufl. Springer-Verlag Berlin Heidelberg New York.

---

1. Der Name \(W\) für die Zustandsvariable ist willkürlich (\(W\) wie Wasserstand). T für Tiefe geht nicht, weil T meistens für die Temperatur benutzt wird, und P für Pegel kollidiert mit dem in aquatischen Ökosystemmodellen üblichen P für Phytoplankton. Derartige Entscheidungen haben oft langanhaltende Auswirkungen.

2. Hier ist \(\left((t-1) \mod 360 \right)\) der Rest bei der Division der Zeit t-1 (in Tagen) durch 360, und \([..]\) symbolisiert das Abrunden auf die nächste ganze Zahl (Gauß-Klammer). Der Tag \(t=400\) liegt also im Februar des zweiten Jahres.

3. Integral von Pegelnull bis Normpegel über die Fläche.