Geoscience Reference
In-Depth Information
they converge towards the same solution, as the spatial differences of k become so
small that every averaging method yields the same result. Harmonic means (not shown
in Figure 9.4 ) underestimate the mean k at the wetting front and the iniltration rate
even more than the geometric mean. In the case of arithmetic averages with ∆ z i = 1 cm,
the calculated iniltration rate is very close to that of the theoretical iniltration curve
(40 compared to 39 mm), while the iniltration rate with geometric averages deviates
more (37 compared to 39 mm). Based on this and other cases, SWAP applies arith-
metic averaging of k and maximum nodal distances at the soil surface of 1 cm.
The second main error source concerned the temporal averaging of water capacity
C . In Question 9.2 the problem is illustrated.
Question 9.2: Suppose at t = j a node in a sandy subsoil has the following state vari-
ables: h i j = -50 cm, θ i j = 0.26210 cm 3 cm -3 and C i j = ∂θ / ∂h = 0.00278 cm -1 . The sub-
soil becomes more dry, and at t = j +1 the node shows the following state variables:
h i j +1 = -52 cm, θ i j +1 = 0.25660 cm 3 cm -3 and C i j +1 = ∂θ / ∂h = 0.00272 cm -1 . The compart-
ment of this node is 5 cm thick. How large is the real difference in water storage of this
compartment? Which water storage difference would you calculate with the expression
C i j ( h i j +1 - h i j ) ∆ z i ? And which with the expression C i j +1 ( h i j +1 - h i j ) ∆ z i ?
Question 9.2 shows the mass balance problem due to temporal averaging of water
capacity C = ∂θ / ∂h . A very elegant solution to this problem was published by Celia
et al. ( 1990 ). Instead of applying during a time step:
( )
θθ
i
j
+
1
−=
j
Ch h
j
+
½
j
+
1
j
(9.3)
i
i
i
i
where C i j denotes some kind of average water capacity during the time step, we may
use the water content estimate at the new time level, θ i j +1, p -1 , in the iterative solution:
(
) +
j
+
1
j
j
+ −
1
,
p
1
j
+
1
,
p
j
+−
11
,
p
j
+−
11
,
p
j
θθ
−=
C
h
h
θ
θ
(9.4)
i
i
i
i
i
i
i
where the superscript p is the iteration level and C i j +1, p -1 is the water capacity evalu-
ated at the pressure head value of the last iteration, h i j +1, p -1 . At convergence, the term
( h i j +1, p - h i j +1, p -1 ) will be small, which eliminates effectively remaining inaccuracies in
the evaluation of C . Implementation of Eq. ( 9.4 ) results in a perfect water balance,
also at larger time steps!
The implicit inite difference solution of Richards' equation that is currently
applied in SWAP therefore reads:
(
) +
j
+− +
11
,
p
j
1
,
p
j
+ −
1
,
p
1
j
+ −
1
,
p
1
j
C
h
h
θ
− =
θ
i
i
i
i
i
j
h
j
+
1
,
p
h
j
+
1
,
p
KK h
j
+
1
,
p
h
j
+
1
,
p
(9.5)
t
z
j
i
1
i
j
j
i
i
+
1
j
j
j
K
+ −
K
t S
i
½
i
½
i
+
½
i
+
½
i
z
z
i
u
Search WWH ::




Custom Search