Evaluation of oil pollution dispersion in an unsaturated sandy soil environment

Document Type : Original Research Paper


1 Ph.D., Graduate in Geotechnical Engineering, University of Kharazmi, Tehran,Iran

2 Associate Professor, University of Tehran, Tehran, Iran

3 Professor, University of Bodenkultur, Wien, Austria

4 Postdoctoral Researcher, University of Brescia, Italy


The current study assesses critical condition of oil dispassion, considering the unsaturated soil condition dispersity behavior for oil dispersion. The numerical model is used as a finite element method to model the oil spill pollution with two different saturated and unsaturated soil conditions chosen and their pollution dispersion results compared. Extracted results from numerical model show that considering the form of unsaturated soil, by changing the matric suction its soil conductivity ratio will differ. Regarding the current study analysis, it has been observed that the pattern of oil dispersion in case of unsaturated soil can be changed, in comparison to saturated soil condition. The vertical penetration of oil pollution in both cases of saturated and unsaturated soil condition will be more than horizontal dispersion pollution speed.  As for oil pollution control in soil domain, the condition of unsaturated soil may be controllable, compared to the saturated one. Extracted results show that oil dispersion velocity, considering saturated soil, is more than 10 times greater than unsaturated one.



An oil spill is the release of a liquid petroleum hydrocarbon into the environment, especially marine areas and soil domain due to human activity, considered a form of pollution (Abbasi Maedeh et al., 2017a). The term is usually applied to marine oil spills, where oil is released into the ocean or coastal waters; however, spills may also occur on land (Abbasi Maedeh et al., 2017b; 2017c). Oil spills may be due to releases of crude oil from tankers, offshore platforms, drilling rigs, and wells, as well as spills of refined petroleum products (such as gasoline, diesel) and their by-products, heavier fuels used by large ships such as bunker fuel, or the spills of any oily refuse or waste oil (Abbasi Maedeh et al., 2016). Accidents that lead to crude oil and refined fuel spills from tanker ships have damaged vulnerable ecosystems in Alaska, the Gulf of Mexico, the Galapagos Islands, France, the Sundarbans, Ogoniland, and many other places. Smaller spills have already proven to have a great impact on ecosystems, such as the Exxon Valdez oil spill due to site's remoteness or the difficulty of an emergency environmental reaction.

Oil spills at sea generally damage as much as those on land. Oil spills on land are more readily containable if a makeshift earth dam can be rapidly bulldozed around the spill site before most of the oil escapes, and land animals can avoid the oil more easily. In general, spilled oil can affect animals and plants in two ways: dirесtly from the oil, itself, and from the response or cleanup process. Animals, using smell to find their babies are unable, due to the strong scent of the oil. This causes a baby to be rejected and abandoned, leaving it to starve and eventually die. Oil can impair a bird's ability to fly, preventing it from foraging or escaping predators. As they preen, birds may consume the oil, coating their feathers, which irritates the digestive tract, alters liver function, and causes kidney damages.

The oil spilling on the earth may occur in two different soil conditions: the saturated or the unsaturated soil condition. An unsaturated soil is usually defined as a three- or four-phase material, consisted of solid, water, air, and the soil-water interface (a contractile skin). Frequently, the last phase is omitted, or considered to be a part of the water phase, due to its negligible thickness (Fig. 1).

Fig. 1. Saturated and unsaturated soil content and condition (Switala, 2016)

Although the presence of a contractile skin is often neglected, it sometimes plays a critical role, especially when the soil is unsaturated significantly. The difference between air pressure and water pressure, acting on the contractile skin, is called matric suction, which is calculated as the following Equation (1) (Switala, 2016).

s = (ua − uw) =      (1)

where ua and uw are air and water pressures, respectively, Ts is the contractile skin surface tension, and Rs is the radius of the contractile skin curvature (meniscus). The Darcy’s law was originally formulated to describe water flow in fully saturated soils (Fig. 2). This velocity is assumed to be proportional to the hydraulic head gradient (one-dimensional flow):

Fig. 2. Natural saturated and unsaturated soil condition and its placement (Switala, 2016)

In unsaturated soils, the coefficient of permeability is often described as a function of saturation degree, or the volumetric water content. In order to calculate the unsaturated coefficient of permeability, kw, the relation of matric suction versus saturation degree is often used (described in the following sections). For either cases the basis for predicting the kw can be the soil pore size distribution (van Genuchten, 1980; Griffiths & Lu, 2005). The coefficient of permeability, obtained from the matric suction against the degree of saturation curve was investigated by Jeff (1999). From the known curve, three soil parameters can be obtained: the air entry value of the soil (ua − uw) b, the residual degree of saturation (Sres), and the pore size distribution index (ps) (Felippa et al., 2001; Fatahi et al., 2010), all of which can be obtained from a graph, shown in Fig. 3, where the saturation is expressed as an effective degree of saturation (Se) (Fatahi, 2007; Fatahi et al., 2009; Buscarnera & Di Prisco, 2012):

In which Sr is the degree of saturation, and Sres is defined as a residual degree of saturation. It is the degree of saturation at which further increasing of matric suction does not make any change in water content of the soil (Zapata et al., 2000). Figure 3 presents the matric suction against the degree of saturation curve.

Fig. 3. Results of matric suction vs degree of saturation, changing in unsaturated soil condition (Switala, 2016)

The following studies have considered the oil dispersion in the soil domain. Lee et al. (2001) studied the effective factors of petroleum contamination propagation in shallow sandy aquifers. Kamon et al. (2004) elaborated the distribution of dense non aqueous phase liquids (DNAPL) in soils, using laboratory and numerical models. Wilson et al. (2006) presented an analytical method to simulate the DNAPLs flow and their transmission in porous mediums. Bandilla et al. (2009) introduced a new method to simulate large-scale subsurface contaminant transport to model reactive contamination transport.

None of the previous studies had investigated the oil dispersion, considering unsaturated soil condition (Jeff, 1999; Geng et al., 2013; Geng & Boufadel, 2015). As such, this study will consider the oil spill dispersion in an unsaturated-saturated soil condition with emphasis on different water table condition and unsaturated soil permeability as well as matric suction. Its finite element theory has been used in previous researches and the theory has been verified by former researchers. The current study's innovation lies in the combination of geotechnical and environmental sciences to assess pollution problems.


Material and Methods

Effective factors in the process of contamination release in the soil include two major processes, consisted of transport processes (also known as advocative- dispersive process) and attenuation ones (Krahn, 2004). Due to the importance of the former, this paper has investigated the main parameters of this process. Two key roles of pollution transportation include pollution shift between two points (advection) and dissipation of pollution during the transport process (dispersion). This research has used general differential equation of pollution propagation and finite element method are in finite element software (Kristine & Vera van, 2014), to make the numerical model and investigate the problem of oil pollution dispersion. The contaminant pollution source is located at the top of the saturated-unsaturated sandy soil surface as well as central line of the model. The characteristic equation of coupled problem is defined as Equation (2) (Mana et al., 2012).


where C is the concentration of oil pollution in a soil domain and specific time; x, a distance variable; t, a time variable; U, a pollution flux speed; D, a hydrodynamic dispersion ratio; θ, the water content of the saturated or unsaturated soil; ρd, the dry density; s, the absorbed pollution; and λ, the half-life ratio of decay. The water content of the entire analysis is assumed to be 0.45 m3/m3 (Reddi et al., 2000). Additional information of numerical model parameters properties are defined as Table 1.

Table 1. Martial properties of numerical model




Soil conductivity


1E-4 to 1E-6

Time of investigation



longitudinal dispersity



Transversal dispersity



Diffusion coefficient



Water table


0 to 5


With regards to the main purpose of the current study, the pollution dispersion would be assessed and the soil absorption and half decay of pollution are neglected. By using the Vanguchten method, the water content of unsaturated soil condition is calculated and reported in Figure 4. 

Fig. 4. Matric suction V.S water content, considering different saturated water content values


Using the unsaturated water content equation and Van Guchten method, the unsaturated soil conductivity considering matric suction and its unsaturated water content is reported in Figure 5.

Fig. 5. Matric suction V.S soil conductivity considering different saturated condition

The 40-meter depth of saturated-unsaturated sandy soil domain is defined as a numerical model boundary, while the infinite boundary is defined for both sides of the numerical model to remove the boxing effects.  The water table of model is changed from 0 to -10 meters. Also, considering the literature, the saturated soil permeability for sandy soil is reported in Table 1. It is assumed that the source of pollution has a point pollution behavior, constructed on top of the soil. The next assumption of this study is constant concentration and flux of oil spilling. The finite element method has been used for this study with the of FEM software to model the spilled oil dispersion (Zou et al., 2013). 


Considering the aim of this study, the mentioned geometry of oil spilling problem is made in numerical software. Generally, the pattern of oil spill with emphasis on time factor around one year is concluded and shown in Figure 6. It has been observed that the pattern of dispersions in both saturated and unsaturated soil is regarded as horizontal and vertical dispersions. To better understand the problem, the first condition of spill assessment concerns saturated soil considering water content of 0.5 and different connectivity of soil, aforementioned. 

Fig. 6. Oil spilling dispersion behavior

With regards to the vertical penetration importance of oil pollution, the depth of penetration of oil (considering the penetration time), and soil connectivity has been plotted in Figures 7 to 9. It has been observed that the water content in saturated soil is not a necessary factor to change the vertical dispersion of oil but the connectivity factor of soil can change the pattern of vertical dispersion of oil. 

Fig. 7. Values of vertical oil concentration, considering K=1e-4 (m/s) in saturated soil condition

Results of Figures 8 and 9 show that the conductivity vaults, being lower than 1e-5 m/s, have had no considerable change in pollution concentration, regarding vertical penetration of oil.

Fig. 8. Values of vertical oil concentration, considering K=1e-5 (m/s) in saturated soil condition

Fig. 9. Values of vertical oil concentration, considering K=1e-6 (m/s) in saturated soil condition

Results of increased water table to 35 m are plotted in Figures 10 to 12. In this case, the top layer of soil is unsaturated while the lower one possesses saturated condition. Concluded results show that the pattern of dispersion in this condition is more different than with saturation soils. It has been observed that the final concentration of oil pollution lately (past one year) is around 10 times lower than saturated soil condition. Since the matric suction can change soil connectivity, the vertical penetration declines to ten times lower that saturated condition. In addition, results show that the dispersion speed of oil spilling during the first days of pollution is lower than saturated soil condition. 

Fig. 10. Values of vertical oil concentration with K=1e-4 (m/s) and water content = 0.5 in unsaturated soil condition

Fig. 11. Values of vertical oil concentration with K=1e-5 (m/s) and water content = 0.5 in unsaturated soil condition

Fig. 12. Values of vertical oil concentration with K=1e-6 (m/s) and water content = 0.5 in unsaturated soil condition

Figure 13 plots horizontal dispersion assessment in case of saturated soil condition. It can be observed that considering 8 m of penetration depth in each horizontal location at earlier times of dispersions, oil concentration near the pollution node is more than 100 times greater than when it is about 10 m far from pollution point center, meaning that the horizontal oil dispersion can be a control, compared to vertical dispersion of oil pollution.

Fig. 13. Values of horizontal oil concentration with K=1e-4(m/s) focusing on saturated soil condition

Figures 14 and 15 give the results of unsaturated soil condition of horizontal dispersion. These results show that the horizontal dispersion in unsaturated soil is determined around 10 times lower that saturated soil in similar condition. 

Fig. 14. Values of horizontal oil concentration with K=1e-4(m/s) focusing on unsaturated soil condition

Fig. 15. Values of horizontal oil concentration with K=1e-6(m/s) focusing on unsaturated soil condition

Considering the assessment results, it has been observed that by dewatering the soil and decreasing the water table condition the horizontal and vertical dispersion can be decreased. Boring, pumping, and draining the water are the most famous techniques to decrease the water table. Additional results show that the speed of dispersion (horizontal and vertical) on the top layers is very faster than the lower ones. By controlling the top layers the pollution can be controlled as well as dewatering technique.


The oil pollution dispersion, in case of saturated and unsaturated soil condition, has been evaluated. Majority of pervious study assumptions involved saturated or dry soil conditions. In case of unsaturated soil condition, the soil conductivity can be changed, considering matric suction of soil condition changing. Results from the current assessment show that the speed of pollution dispersion in saturated soil is far too greater than unsaturated soil condition. The final concentration of oil pollution in case of saturated soil is around 100 times greater than unsaturated soil condition in case of vertical penetration.

The horizontal oil spill dispassion concentration is lower than vertical penetration in both cases of saturated and unsaturated soil conditions. Also results show that to decrease the speed of horizontal and vertical dispassion changing, the soil condition from saturate to unsaturated soil is the most effective method. The rate of penetration and dispersion of oil pollution during the first two months of saturated and unsaturated soil are more different, but by increasing the dispassion time their behavior would be the same in different concentrations.


The authors wish to appreciate the institute of geotechnical engineering, University of Bodenkultur, Wien, Austria, for their generous supports to complete this work.

Abbasi Maedeh, P., Ghanbari, A. and Wu, W. (2017a). A new analytical model for natural period analysis of elevated tanks considering fluid-structure-soil interaction. J. of Geoeng., 12(1): 1-12.
Abbasi Maedeh, P., Ghanbari, A. and Wu, W. (2017b). Estimation of elevated tanks natural period considering fluid-structure-soil interaction by using new approaches. Earthq. Struct., 12(2): 145-152.
Abbasi Maedeh, P., Ghanbari, A. and Wu, W. (2017c). Investigation of soil structure interaction and wall flexibility effects on natural sloshing frequency of vessels. Civil Eng. J., 3(1): 45-56.
Abbasi Maedeh, P., Ghanbari, A. and Wu, W. (2016). Analytical assessment of elevated tank natural period considering soil effects. Coupled Systems Mechanics, 5(3): 223-234.
Bandilla, K.W., Rabideau, A.J. and Jankovic, I. (2009). A parallel mesh-free contaminant transport model based on the Analytic Element and Streamline Methods. J. Adv .Water Resour., 32: 1143-1153.
Buscarnera, G., and Di Prisco, C. (2012). Discussing the definition of the second-order work for unsaturated soils. Int. J. Numer. Anal. Meth. Geomech., 36: 36-49.
Fatahi, B. (2007). Modelling of influence of matric suction induced by native vegetation on sub-soil improvement. PhD thesis, University of Wollongong, Australia.
Fatahi, B., Khabbaz, H. and Indraratna, B. (2010). Bioengineering ground improvement considering root water uptake model. Ecol. Eng., 36: 222-229.
Fatahi, B., Khabbaz, H. and Indraratna, B. (2009). Parametric studies on bioengineering effects of tree-root based suction on groud behaviour. Ecol. Eng., 35: 1415-1426.
Felippa, C.A., Park, K.C., and Farhat, C. (2001). Partitioned analysis of coupled mechanical systems. Compt. Meth. Appl. Mech. Eng., 190: 3247-3270.
Geng, X. and Boufadel, M.C. (2015). Numerical study of solute transport in shallow beach aquifers subjected to waves and tides. J. Geophys. Res. Oceans, 120(2): 1409-1428.
Geng, X., Boufadel, M.C. and Wrenn, B. (2013). Mathematical modeling of the biodegradation of residual hydrocarbon in a variably-saturated sand column. Biodegradation, 24(2): 153-163.
Griffiths D.V. and Lu, N. (2005). Unsaturated slope stability analysis with steady infiltration or evaporation using elasto-plastic finite elements. Int. J. Num. Anal. Meth. Geomech., 29: 249-267.
Jeff, K. (1999). Practical design calculations for groundwater and soil remediation. Boca Raton, CRC Press LLC: 152-163.
Kamon, M.K., Junichi, I.T. and Katsumi, T. (2004). Two-dimensional DNAPL migration affected by groundwater flow in unconfined aquifer. J. Haz. Mat., 110: 1-12
Krahn, J. (2004). Seepage modeling with SEEP/W: An engineering methodology, GEO-SLOPE International Ltd., Calgary, Alta., Canada.
Kristine, V. and Vera van, B. (2014). 3D finite element method (FEM) simulation of groundwater flow during backward erosion piping, Front. Struct. Civ. Eng., 8: 160-166.
Lee, J.Y., Cheon, J.Y., Lee, K.K., Lee, S.Y., and Lee, M.H. (2001). Factors affecting the distribution of hydrocarbon contaminants and hydro geochemical parameters in a shallow sand aquifer. J. Cont. Hydrol., 50: 139-158.
Mana, D.S.K., Gourvenec, S.M., Randolph, M.F. and Hossain, M.S. (2012). Effect of gapping on the uplift resistance of a shallow skirted foundation. In: Proc. 7th Int. Offshore Site Inv. & Geotech. Conf., Society of Underwater Technology, London: 403-410.
Reddi, L.N., Lee, I.M. and Bonala, M.V.S. (2000). Comparison of internal and surface erosion using flow pump tests on a sand-kaolinite mixture. Geotech. Test. J. 23: 116-122.
Switala, B.M. (2016). Analysis of slope stabilization by soil bioengineering methods. PhD thesis, University of Natural Resources and Life Sciences, Vienna.
van Genuchten, M.T. (1980). Academia form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J., 44(5): 892-898.
Wilson, C.S., Weaver, J.W. and Charbeneau, R.J. (2006). A screening model for simulating DNAPL flow and transport in porous media: theoretical development. Environ. Model. Soft., 21: 16-32.
Zapata, C.E., Houston, W.N., Houston, S.L., and Walsh, K.D. (2000). Soil water characteristic curve variability. In C.D Shackelford, S.L. Houston, and N.Y. Chang, editors, Advances in Unsaturated Geotechnics, Denver, USA, Geo-Institute. 118.
Zou, T.H., Chen, Q., Chen, X.Q. and Cui, P. (2013). Discrete numerical modeling of particle transport in granular filters. Comput. Geotech., 47: 48-56.