Das Modell wurde von ERROR: Zotero entry with cite key Allen2003 not found. beschrieben. Das Xinanjiang-Modell beruht auf der grunsätzlichen Idee, den im Einzugsgebiet vorhandenen Speicher und den Anteil der gesättigten Flächen in Abhängigkeit der vorhandenen Feuchte zu beschreiben.

Der Kern des Modells ist eine Funktion zwischen der gesättigten Fläche $A_s$ und der aktuell im Oberboden gespeicherten Feuchte $i$. Für jede im Oberboden vorhandene Feuchte kann eine äquivalente Grenze $i_o$ für die Sättigung und eine entsprechende gesättigte Fläche $A_s$ angegeben werden. Für das Gebiet gibt es eine maximal mögliche \index{Speichermenge} Speichermenge für Feuchte im Oberboden von $i_m$.

$$i = i_m*[1-(1-A_i)^{1/b}]$$

Die Funktion könnte auch über eine hypsometrische Kurve oder über eine kumulierte Verteilungsfunktion der nutzbaren Feldkapazität (Verdunstung) bzw. der Feldkapazität (Sickerung) angenähert werden. Die Versiegelung könnte über einen festen Flächenanteil eingebaut werden, ebenso Seeflächen.

In R programmiert ergibt sich für den Bereich $0 < A < 1$ und einen maximalen Speicher $i_m$ folgende Funktion.

VIC.R

# Flaeche vs. Saettigung
A <- seq(0.0,1,0.05) # Flaeche zwischen 0 und 1
d <- A*0.0           # Defizit initialisiert
im <- 1.0            # maximale Speicherung, 0 bis 1
b <- 0.5             # Formfaktor
io = 0.5             # aktuelle Feuchte
i <- im*(1-(1-A)^{1/b})
for (n in 1:(length(A))) {
    if (i[n]>io) {
        d[n] <- (i[n]-io)
    } else {
        d[n] <- 0.0
    }
}
plot(A,i)
curve((1-(1-x)^(1/1.5)),0,1,
    add = TRUE, col = "blue")
curve((1-(1-x)^(1/1)),0,1,
    add = TRUE, col = "red")
curve((1-(1-x)^(1/b)),0,1,
    add = TRUE, col = "gray")
curve((1-(1-x)^(1/0.2)),0,1,
    add = TRUE, col = "green")

Für eine konkrete Speichermenge $i_o$ kann nun das jeweils pro Fläche vorhandene Defizit und die gesättigte Fläche berechnet werden.

$$d=i-i_o = i_m*[1-(1-A_i)^{1/b}]-i_0$$

Ae.R

# Flaeche vs. Saettigung
A <- seq(0.0,1,0.05) # Flaeche zwischen 0 und 1
d <- A*0.0           # Defizit initialisiert
im <- 1.0            # maximale Speicherung, 0 bis 1
b <- 0.5             # Formfaktor
io = 0.5             # aktuelle Feuchte
i <- im*(1-(1-A)^{1/b})
for (n in 1:(length(A))) {
    if (i[n]>io) {
        d[n] <- (i[n]-io)
    } else {
        d[n] <- 0.0
    }
}
plot(A,d,xlim=c(0:1),ylim=(0:1))
curve((1-(1-x)^(1/1.5))-io,0,1,
        add = TRUE, col = "blue")
curve((1-(1-x)^(1/1))-io,0,1,
        add = TRUE, col = "red")
curve((1-(1-x)^(1/b))-io,0,1,
        add = TRUE, col = "gray")
curve((1-(1-x)^(1/0.2))-io,0,1,
        add = TRUE, col = "green")

Das übrige hydrologische Verhalten ergibt sich aus dieser Speicherfunktion. Der Niederschlag, der auf die gesättigte Fläche fällt, wird sofort Abfluss. Der Niederschlag, der auf die ungesättigte Fläche fällt, infiltriert mit einer hydraulischen Leitfähigkeit, die sich aus dem Defizit über die Brooks-Corey oder die Van-Genuchten-Beziehung ergibt. Die entsprechenden Speicheränderungen werden berücksichtigt und der Speicherinhalt wird nach dem Zeitschritt aktualisiert. Der Speicherinhalt wird über die Wasserhaushaltsgleichung aktualisiert. Dabei ist $i_{i+1} = i + F - E- S$

Über den Speicherzustandd $i$ wird die gesättigte Fläche aktualisiert. Dafür gilt folgende Formel, die sich Umstellung der Funktion $A(i)$ ergibt.

Auch die Drainage/Kanalisation oder Makroporen können/sollten berücksichtigt werden, indem die entsprechenden Flächen eine schnelle Abflussbildung oder einen schnellen Transfer in das Gewässer unter Umgehung des Einzellinearspeichers enthalten. Der ELS könnte durch eine Boussinesq-Funktion ersetzt werden oder durch eine Speicherkaskade.

$$ A_i = 1 - \left(1-\frac{i}{i_m}\right)^{b} $$

S.R

# Flaeche aus Saettigung
im<- 1.0            # maximale Speicherung, 0 bis 1
i <- seq(0.0,im,0.05)# Flaeche zwischen 0 und 1
b <- 0.5             # Formfaktor
io<- 0.5             # aktuelle Feuchte
A <- 1 - (1-(i/im))^b
plot(i,A,xlim=c(0:im),ylim=(0:1))
curve(1-(1-(x/im))^b,0,im,
        add = TRUE, col = "gray")
curve(1-(1-(x/im))¹,0,im,
        add = TRUE, col = "red")
curve(1-(1-(x/im))¹.5,0,im,
        add = TRUE, col = "blue")

Die Verdunstung von der gesättigten Fläche entspricht der potentiellen Verdunstung. Somit wird implizit auch kapillarer Aufstieg berücksichtigt. Die Verdunstung aus der ungesättigten Fläche wird über die Brooks-Corey oder die Van-Genuchten-Beziehung reduziert. Hier ist die Berechnung nach der Brooks-Corey dargestellt.

$$E_a = E_{pot} * [(i-i_r)/(i_{max}-i_r)]^{λ}$$

Die Sickerung aus dem Oberboden in das nächst tiefere Kompartment richtet sich ebenfalls nach der Sättigung und wird über die Brooks-Corey oder Van-Genuchten-Beziehung berechnet. Hierbei kann eine vollständige Richards-Beziehung mit

$$K(θ) = K_s * [(i-i_r)/(i_{max}-i_r)]^{γ}$$

Aus dem Wert $K(θ)$ wird mit der Richards-Gleichung die Sickerung berechnet - dafür muss der Gradient zwischen dem Oberboden und dem unteren Bodenspeicher/Grundwasserkompartment berücksichtigt werden.

$$v_f = k_f(θ(ψ))*\frac{d(h+ψ)}{dz}$$

Die entscheidende Annahme ist, dass es eine funktionale Beziehung zwischen dem Defizit $d$ und der Bodenfeuchte $\theta$ bzw. der Saugspannung $\psi$ gibt, die durch eine pf-Kurve oder Sättigungs-Saugspannungsbeziehung beschrieben werden kann. Diese wird mit einer Brooks-Corey oder einer Van-Genuchten-Beziehung dargestellt und ist Grundlage der Anpassung von $E_a/E_{pot}$ und von $v=-k_f(θ(ψ))$.

Der Abfluss aus dem tieferen oder Grundwasserspeicher erfolgt mit einem Einzellinearspeicher. Ein Einzellinearspeicher leert sich nach der mathematischen Beziehung:

$$h_i = (h_{i-1}-h_B)*\kappa$$

Hierbei ist $\kappa$ eine Rezessionskonstante. $h_{i-1}$ ist der Wasserstand des vorigen Zeitschrittes und $h_B$ ist der Wasserstand des Vorfluters, der nicht unterschritten werden kann. Diese differenzielle Form wird bei kleinen Zeitschritten verwendet. Für längere Zeiträume gibt es eine analytische Lösung:

$$h = h_0*e^{-\kappa*t}$$

Hierbei ist $\kappa$ wieder die Rezessionskonstante, $h_0$ ist der Anfangswasserstand von dem ein exponentieller Abfall über die Zeit $t$ erfolgt.

Die Umwandlung des Oberflächenabflusses als Einzugsgebietsantwort erfolgt mit einem Unit Hydrograph. Hier kann der geomorphologic Unit Hydrograph verwendet werden.

Die Weiterleitung des Abflusses zum einem flussabwärts gelegenen Pegel erfolgt mit dem Muskinghum-Cunge-Routing-Verfahren.

Die Verdunstung kann mit Penman-Monteith berechnet werden, wenn es nicht genügend Daten gibt, kann Priesley-Taylor, Hamon oder eine einfachere Verdunstungsformel verwendet werden.

Nach Bedarf können Seen in das Modell eingebaut werden.

Die Schneeschmelze könnte mit einem sehr ähnlichen Ansatz berücksichtigt werden oder mit einem einfachen Gradtag-Verfahren. Durch Diskretisierung können unterschiedliche Bodenarten berücksichtigt werden, indem Teileinzugsgebiete mit homogenem Boden berücksichtigt werden. Die Landnutzung beeinflusst die für den Transport angesetzte Bodentiefe, die potentielle Verdunstung, den Interzeptionsspeicher und Schneeschmelze.

| victhl.r

# Parameter für V
A = 2
D = 180
S = 2.5
 
# Verdunst. erz.
V = function(A=2,D=90,S=2.5,Ta){
  E = A*sin((Ta+D)/365*2*pi)+S
  return(E)
}
 
Ta = seq(1,365,1)
N = numeric(length(Ta))
Q = numeric(length(Ta))
I = numeric(length(Ta))
A = numeric(length(Ta))
Ea = numeric(length(Ta))
 
# Parameter
imax = 100
imin = 5
ih = 10
i0 = 20
A0 = 5
lambda = 0.5
eta = 2
beta = 1.5
 
# initialisieren
Ah = A0
ih = i0
n = 0
 
for (ti in Ta) {
 
  n = n+1
  N[n] = runif(1,0.1,2)
  Ep=V(A=2,D=90,S=2.5,ti) 
 
  # Defizit
  if(n==1){
    A[n] = Ah
    I[n] = ih
    Q[n] = 0
    Ea[n] = Ep
    def = imax-ih
  } else{
    def <- imax*(1-(1-A[n])^{1/beta})-I[n]
    # Verdunst rechnen
    Ea[n] = Ep*((I[n]-imin)/(imax-imin))^eta
  }
 
  # Änd. der Feuchte
  di = N[n]-Ea[n]
  if (I[n]+di<def){
    I[n] = I[n]+di
  } 
  else{
    I[n] = imax
    Q[n] = imax-(I[n]+di)
  }
 
  # neuer Flächenanteil
  A[n] <- 1 - (1-(I[n]/imax))^beta
 
  # Abfluss rechnen
  Q[n] = N[n]*A[n]
}
hydro = data.frame(N,Ea,A,I,Q)