Development of a numerical solution for transport simulation in unsaturated porous media

Groundwater contamination has become increasingly prevalent in recent years, emerging as a major environmental issue. Effective groundwater management requires tools that can track and predict the movement and behaviour of contaminants within subsurface environments. There is an increasing reliance on using numerical groundwater modeling to predict fate and transport of contaminants as well as employing predictive simulations to anticipate concentrations of contaminants in groundwater in order to evaluate a remedial action. Soils are porous media and exhibit a heterogeneous, multi-scale architecture. Their non-linear material properties have a significant influence on the soil water dynamics, which poses difficulties for numerical models. This study examines the numerical approximations applied to model water flow within the variably saturated soil zone, and transport of particles through this unsaturated condition. The governing equation that describes the flow of water in unsaturated soil profile is Richards’ equation which is a degenerate parabolic nonlinear partial differential equation. Getting robust, accurate and cost-effective results is difficult in particular because of moving sharp wetting fronts due to the nonlinear hydraulic properties. A numerical model, using the finite difference method for formulating the equation is used to simulate water flow in unsaturated soil. The approach is checked through various test-cases that incorporate both constant head, constant flux, and free drainage boundary conditions, along with soils exhibiting various functional forms for moisture content and hydraulic properties. Then, the result of simulation is feed to a validated solute transport model developed at Politecnico di Torino to observe the transport of particles in the vadose zone. Numerical simulations are compared with the results of Hydrus1D software to assess the effectiveness of the modeling and validation of the proposed code. To this aim, statistical measures (NSE, and R2) are employed and the modeling results illustrate acceptable values. The Code has been written in MATLAB R2023b platform and run in an Intel Corei7-1075H CPU @2.60GHz system.