A Time-Splitting technique for the solution of Density Dependent Flow and Transport in Groundwater

The mathematical model of density-dependent flow in groundwater can be formulated as a coupled system of two partial differential equations, one describing mass conservation for the water-salt solution (the flow equation), and the other mass conservation for the salt contaminant (the transport equation) (Gambolati, Putti, & Paniconi 1999). From a physical point of view, the salt contained in the groundwater affects the solution density and induces changes in the flow field, which becomes dependent not only on the hydrogeological parameters of the aquifer, but also on the salt concentration. This interdependence produces a coupling between the flow and transport equations. The solution of the discretized system of equations is usually addressed with an iterative Picardlike scheme by which the problem is decoupled by first solving the flow equation, then calculating the velocity field, and finally solving the transport equation. This three-step sequence is repeated until convergence. A difficulty for the successful application of the above procedure is the requirement that accurate velocity fields be obtained from the solution of the flow equation. Another problem that may influence the convergence of the Picard method is oscillatory behavior of the solution when the transport equation is advection dominated (Huyakorn, Andersen, Mercer, & White 1987), (Putti & Paniconi 1995). To overcome the previous problems, an algorithm is developed, based on the Mixed Hybrid Finite Element (MHFE) method for the flow equation (Bergamaschi & Putti 1999), and a combination of MHFE with high resolution Finite Volume (FV) scheme for the transport equation. In the latter approach, discretization in two spatial dimensions is obtained by means of triangle-based lowest order MHFE and TVD (Total Variation Diminishing) FV for the dispersive and advective fluxes, respectively. MHFE and FV are known to be stable and accurate (both are globally second order accurate in space). Combination of these two schemes is implemented via a time-splitting algorithm (Mazzia, Bergamaschi, & Putti 2000; Mazzia, Bergamaschi, & Putti 1999) allowing an explicit timestepping for FV and an implicit time-stepping for MHFE. This algorithm is a generalization to unstructured grids of the Godunov-Mixed Methods (GMM) developed by (Dawson 1993; Dawson 1995). The overall procedure is shown to attain first order accuracy in time, even though a second order Runge-Kutta time stepping is used. To resort to global second order accuracy a correction term is added in the FV equation (Mazzia, Bergamaschi, Putti, & Dawson 2000).