Notes: This section, as an excerpt, refers to figures not presented on this web site. You must order the full publication to see them. It will be helpful to go to Gathered Equations and print out the equations before reading this section. This page shows greek characters correctly with Netscape Navigator 3.0 on a Win 95 computer with the SYMBOL TrueType font set. If it does not do so on your computer, then psi may show as y, eta as h, and delta as D.
Just enough material has been presented here to make the point that Darcian imbalance errors in some numerical methods can predict backwards flows. More material, plots and examples are available in the full report.
The reader should be aware that the errors presented thus far are not the only ones of concern. It is common enough for unsaturated zone modelers to separate flow into the diffusive flow controlled by Kh and the gravity flow controlled by Kx. Consider the implications to equation [6]. If Kh and Kx have errors or are simply not Darcian, then there is a singularity in [6] about Dx à Dy that can drive the resulting Kv to large positive or negative values. Since these are nearhydrostatic conditions, it is reasonable to expect this effect in modeling a perched water table or horizontal capillary barrier, say, a horizontal fracture in Topopah Spring welded tuff, where the vertical downwards flow is backing up.
This makes it reasonable to ask: How many modelers make a serious effort to assure mass balance to one part in 10^{6} without balancing equation [6]? What is the magnitude of the error? How can it be corrected?
Zaidel and Russo (1992) developed a mean for diffusive flow that was very much like [45], except that the premise, if one remembers correctly, was the integral of a power function of (a+bx) in [7]. But for Kx they used the arithmetic mean, ka [18]. For the moment, let us follow their example as one which modelers may reasonably be using with BrooksCorey relations. But let us use the correct, Darcian Kh in [6], which is readily derivable, and substitute some nonDarcian means for Kx to get [52]. Let us use two standard means for k*, such as ka and kg [19], as well as the two that we have identified here as the limits of transformation with u and Dx, kxo = km(m=11/h) [45] and kr2. Let the resulting estimates of Kv, divided by Kv, be called kza, kzg, kzxo and kzk2.
[52]
Figures 32ad use the data for Kv in Figures 25ad to plot the ratio kza in the plane of (kr2/kr1) and Dx/Dyd, or "scaleddx", for h = 2/3, 1.9825, 8.075 and 20. Similarly, Figures 33ad demonstrate kzg, 34ad demonstrate kzxo, and 35ad demonstrate kzk2. Here we have converted the scaling parameter, u, to Dx/Dyd via [49] to make the plots more accessible to modelers who cannot be expected to have any familiarity with u.
































Note that most of these figures are plotted at linear vertical scales that show only the gross errors. Where points seem to be missing from the plot, they are either off the range positive or negative. Table 6 shows the maximum and minimum values of the plotted data. It has no inherent value except to illustrate the possible magnitudes of the Kvestimate errors near hydrostatic conditions.
Figure 33b is plotted on a log vertical scale because it has no negative values. Many of the plots seem to be ripped apart as a strong man might tear a phone book, demonstrating the discontinuity in [45] and [52], for dryoverwet conditions and large space steps. The region of tearing is the division between upwards and downwards flow, at hydrostatic conditions.
Added material: In the regions where the estimates of Kv are negative, those estimates predict backwards flow. The full report shows examples for a composite conductivity relation, derived from the assumption that matrix and fracture flow are both vertical and parallel in a supposed welded volcanic tuff. These examples show how even the estimate of the hydrostatic conditions can be misplaced in the (kr1,kr2)plane by Darcian imbalance errors. In those cases, the flow can be backwards either upwards or down, depending on errors and conditions. This has implications any such modeling about hydrostatic conditions, from capillary barriers to perched water tables to the transition between infiltration and evaporation. DLB




































































Are the Darcian mean approximations presented here for the BrooksCorey relation subject to this kind of error? Yes, and unfortunately this investigator is not yet prepared to offer up a fix, other than to suggest some possible approaches.
..... material deleted .....
One other possibility is to use some other numerical method that does not involve equation [6]. For example, it is possible to formulate a mixed variable finite difference solution of Richards' equation [57] using water content, q, in the time term and total head, H = xy, in the space term [58]. This method is at least superficially diffusive, without explicit reference to the gravity flow, allowing the use of only Kv (for K_{i+1/2} & K_{i1/2}). The reader may find, however, that developing a Jacobian matrix for Newton's method, using [58] and approximations of Darcian means, to be a challenge in stability. And, of course, the approximations for Kv are a bit worse than for Kx.
[57]
[58]
Nevertheless, if this investigator is right in his assessment of the importance of Darcian means, then, at the least, no exposition of an interblock conductivity mean for vertical unsaturated flow, given a particular conductivity relation, is complete without an analysis of the errors similar to that presented here. Mass balance is insufficient to guarantee flow accuracy unless the interblock conductivities can be demonstrated to be both reasonably Darcian and properly balanced among the total, diffusive and gravity flows.
For who would want to stand up in court to defend an environmental impact statement, only to have opposing counsel point out that some of the calculated flows are backwards? It could well become incumbent for those who would use nonDarcian means to demonstrate in every case that they are indeed benign. For if one points out that a method of adding 1+1 does not produce 2, then those who use the method have the responsibility of demonstrating the validity of all the cases.
REFERENCE:
Zaidel, D. and D. Russo, 1992, Estimation of finite difference interblock conductivities for simulation of infiltration into initially dry soils, Water Resour. Res. 28(9):22852295.