Two Dimensional Numerical Modelling Of Variably Saturated Flows ?

Abstract

The prediction of moisture and contaminant transport through unsaturated soil to groundwater is becoming increasingly important in the fields of hydrology, agriculture, and environmental engineering. Computer-aided simulation techniques enable one to conduct a series of systematic numerical experiments to analyze flow phenomena in subsurface hydrology under various physical and chemical processes. The flow movement depends upon medium characteristics, initial and boundary conditions, which reflect physical processes occurring below the ground. To understand the effects of physical processes, an efficient and accurate model is needed. Thus, the model developed must be able to handle varied initial and boundary conditions. In this regard, infiltration into very dry soil becomes a very important problem of study.

Most of the earlier numerical models developed are concentrated on the development of an efficient algorithm or the modeling of a particular process which governs the flow in unsaturated or saturated-unsaturated homogeneous medium. Not much work has been done on the analysis of variably saturated flow in layered soil media. Models to simulate unsaturated flow through dry soils, especially through layered soils with varied boundary conditions, are very limited. Further, not many studies have been reported in the literature on the prediction of seepage face development and the phreatic surface movement in variably saturated media with layering. These aspects are very important in determining the flow field and the discharge from the domain. A detailed literature review covering the above aspects has been made and is reported in this thesis.

In the present work, two-dimensional numerical models to predict the movement of the wetting front in the unsaturated domain and the movement of the phreatic surface in homogeneous and layered porous media under various initial and boundary conditions are developed based on finite difference and finite volume techniques. These models can handle flow in both rectangular flow domains and radial flow domains. The initial condition settings include the handling of very dry soil medium without any transformation of the governing equation, handling of infiltration and constant head conditions at the boundaries under steady state as well as transient scenarios. The models are also able to handle various soil moisture characteristics which depict the nonlinear behavior between hydraulic conductivity, moisture content, and pressure head in a soil medium.

A mixed form of the governing partial differential equation is used in the present study as it leads to better mass conservation. The finite difference model uses a central difference approximation for the space derivatives and an Eulerian backward difference approximation for the time derivative. The fully implicit formulation is solved iteratively using the Strongly Implicit Procedure after making the Picard approximation for the highly nonlinear coefficients. The process of infiltration into an initially dry soil leads to the development of a steep wetting front. As the finite volume technique is naturally an upwind method, it is expected to play a positive role in modeling such processes accurately. Hence, a finite volume model is also developed by approximating the convective part using a MUSCL approach and a fully implicit central difference method for the diffusive part of the governing equation.

The models developed are validated using both experimental data and numerical solutions for problems reported in the literature. The validation problems cover a wide range of physical scenarios such as: infiltration into a very dry soil, infiltration into a dry soil column with gravity drainage, development of water table mound, steady state drainage in a sand-filled wedge-shaped tank with seepage face development, and transient seepage face development in a rectangular domain. Five test problems are used for the validation of the models. The developed models perform very well for the test problems considered, indicating the models' capability in handling such situations. The results obtained by using the present models for simulating flow through highly unsaturated (very dry) soils show that the models perform very well when compared with models that use transformation techniques to handle such problems. The performance of the present models in comparison with the experimental data and numerical models available in the literature shows the suitability of the present models in handling such situations.

The present models are also used to analyze various types of unsaturated flow problems with varying initial and boundary conditions. The boundary conditions considered are no flow and/or free flow conditions along the left, right, and bottom boundaries with infiltration condition along a part of the top boundary. For the various cases considered in the present study, the infiltration rate varies from 5 cm/day to 50 cm/day through an initially very dry soil of —15000 cm pressure head in homogeneous and layered soils. Different types of soil media considered vary from sandy loam, loam, and clay with horizontal and vertical layering of these soils. A total number of 14 cases are analyzed. The results are discussed in terms of pressure head variation in the flow domain along with moisture redistribution for all the cases under consideration. It is observed from these studies that the infiltration rate plays an important role in the wetting front movement through layered soils depending on the type of layering and the boundary conditions considered. The soil properties of various layers affect the movement of the wetting front by changing the direction of movement. Even though the wetting front movement is predominantly vertical, there is a tendency for the wetting front to move in the horizontal direction as it moves from a coarse soil to a fine soil. It is also observed that the vertical layering of soils with different hydraulic conductivity helps in redirecting the flow towards the bottom boundary through the neighboring coarser layers.

As the finite volume method is more suitable for simulating sharp fronts, it is expected to perform better than the finite difference method for simulating infiltration into very dry soils. So, a comparison is made between the performance of these two models by using the above test problems. It is observed from these studies that the performance of both models is the same except that the finite volume method takes much more CPU time than the finite difference model. Considering the type of problems tested, it is observed that for modeling unsaturated and saturated-unsaturated flows, the finite difference method is better in comparison to the finite volume method. It may be due to the predominant diffusive characteristics of the governing equation even while modeling flow through very dry soils.

Proper estimation of the seepage height is an important aspect in finding the discharge through the porous medium. It is observed from the literature that the use of a saturated flow model in such situations can lead to an underestimation of the discharge through the porous medium. This effect is more important when dealing with small dimension problems. It is also observed that various parameters which govern the moisture movement through saturated-unsaturated regions affect proper estimation of the seepage face height and thereby discharge. Various factors like the effect of boundary conditions, type of soil layering, problem dimension, and aspect ratio on seepage face development and the associated phreatic surface formation are studied in the present work. It is seen from the present study that the seepage face development is more in rectangular flow domains than in radial flow domains for both homogeneous and layered soils. It is also seen that the seepage face development in rectangular problems is more sensitive than radial flow problems for various factors considered. The seepage height is also influenced by the tailwater level. It is seen from the present study that as the tailwater level increases, the seepage face reduces with no seepage face development for some of the cases studied. This influence is relatively less for radial flow problems. As the length of the domain increases, the seepage height decreases. It is seen that for different cases with the same horizontal dimension, as the height of the domain increases, the seepage face height also increases. This phenomenon is observed for both homogeneous and layered soil mediums. The influence of the aspect ratio, which is the ratio of the length to height of the domain, indicates that as the aspect ratio increases, the seepage height decreases. The type of the soil layering is observed to have a very strong influence on the seepage face development. The study for understanding the effect of soil layering on the development of the seepage face and phreatic surface suggests that as the coarseness of the material increases, the phreatic surface becomes flatter, and its steepness increases with the fineness of the soil.

The present model is also used for studying the transient phreatic surface movement and the seepage face development. This is studied for homogeneous and layered rectangular soil media. The present study is used to understand the effect of specific storage on the phreatic surface movement and the seepage face development. The studies indicate that the influence of specific storage on the seepage face development is insignificant in homogeneous soils with only very little effect in the early time for longer domains. It is also observed that the influence of specific storage is significant in the case of layered soils. This effect depends on the type of layering and the problem dimension and is observed to have influence for a relatively longer period. This observation suggests the importance of specific storage on transient seepage face development. When the specific storage effect is considered, the drainage of the soil becomes faster, resulting in a faster decline of the phreatic surface with time. The influence of specific storage is also studied considering the problem dimension effect. It is seen that as the aspect ratio increases, the effect of specific storage on the phreatic surface development decreases. The studies with a change in the upstream boundary condition from a constant head to a no flow condition indicate that the effect of specific storage has no significant influence on the phreatic surface development for both homogeneous and layered soils.