Report DOE/ER/82329-2, Backwards Flow Errors

prev   index   next

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 near-hydrostatic 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 106 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 Brooks-Corey relations. But let us use the correct, Darcian Kh in [6], which is readily derivable, and substitute some non-Darcian 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=1-1/h) [45] and kr2. Let the resulting estimates of Kv, divided by Kv, be called kza, kzg, kzxo and kzk2.




Figures 32a-d use the data for Kv in Figures 25a-d to plot the ratio kza in the plane of (kr2/kr1) and Dx/Dyd, or "scaled-dx", for h = 2/3, 1.9825, 8.075 and 20. Similarly, Figures 33a-d demonstrate kzg, 34a-d demonstrate kzxo, and 35a-d 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.



 Figure 32a: kra for h = 2/3


 Figure 32b: kra for h = 1.9825


 Figure 32c: kra for h = 8.075


 Figure 32d: kra for h = 20



 Figure 33a: krg for h = 2/3


 Figure 33b: krg for h = 1.9825


 Figure 33c: krg for h = 8.075


 Figure 33d: krg for h = 20



 Figure 34a: kzxo for h = 2/3


 Figure 34b: kzxo for h = 1.9825


 Figure 34c: kzxo for h = 8.075


 Figure 34d: kzxo for h = 20



 Figure 35a: kzk2 for h = 2/3


 Figure 35b: kzk2 for h = 1.9825


 Figure 35c: kzk2 for h = 8.075


 Figure 35d: kzk2 for h = 20


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 Kv-estimate 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 dry-over-wet 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


Table 6: Minimum and Maximum Values of Kv-estimates for Figs. 32-35






























































Are the Darcian mean approximations presented here for the Brooks-Corey 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 = x-y, 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 Ki+1/2 & Ki-1/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.





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 non-Darcian 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.



 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):2285-2295.


This page was designed by Aquarius Engineering © 1998