Moisture Transfer Finite Element Modeling with Soil-Water Characteristic Curve-Based Parameters and its Application to Nhan Co Red Mud Basin Slope

VNU Journal of Science: Earth and Environmental Sciences, Vol. 37, No. 1 (2021) 103-115 
103 
Original Article 
Moisture Transfer Finite Element Modeling with Soil-Water 
Characteristic Curve-Based Parameters and its Application to 
Nhan Co Red Mud Basin Slope 
Nguyen Van Hoang1,*, Hoang Viet Hung2, Pham Van Dung2 
 1Institute of Geological Sciences-Vietnam Academy of Science and Technology, 
84 Chua Lang street, Lang Thuong, Dong Da, Hanoi, Vietnam 
2Thuyloi University - Hanoi Campus, 175 Tay Son, Dong Da, Hanoi, Vietnam 
Received 16 July 2020 
Revised 14 December 2020; Accepted 13 February 2021 
Abstract: Since the year of 2017, landslides at the red mud basins in Nhan Co alumina factory, Dak 
Nong province have been occurring during the rainy season. The change of the soil physical and 
mechanical parameters due to rainwater infiltration has been considered as the main factor of the 
slope instability. The soil cohesion and angle of internal friction depend greatly on the soil moisture. 
Specifically, soil with a lower moisture content has a higher shearing strength than that in soil with 
higher moisture content. The finite element modeling of moisture transfers in unsaturated soils 
through the relationship between soil moisture, soil suction, unsaturated permeability and soil-
moisture dispersivity is capable of accurately predicting the wetting front development. The element 
sizes and time steps have been selected based on detailed analysis of analytical error estimation and 
on the numerical simulations with different element sizes numerical simulation errors. Soil samples 
had been taken then the soil different suctions and corresponding soil moisture values have been 
determined in the laboratory. The soil water characteristic curve (SWCC) parameters (a, n and m) 
have been determined by the best fitting using the least squared error method. The hydraulic 
conductivity of the saturated soil, one of the key input parameters was also determined. The results 
of the application to the study area's slope has shown that the wetting front depth could be up to 8 
meters for 90 days of moisture transfer due to the rainwater infiltration. The wetting front depth and 
the length of the intermediate part of the moisture distribution curve have increased with the 
infiltration time. The soil moisture distribution with a depth is an essential information to have soil 
strength parameters for the slope stability analyses. The slope stability analysis with the soil shear 
________ 
* Corresponding author. 
 E-mail address: VDC@yahoo.com 
 https://doi.org/10.25073/2588-1094/vnuees.4655 
N.V. Hoang et al. / VNU Journal of Science: Earth and Environmental Sciences, Vol. 37, No. 1 (2021) 103-115 104 
strength parameters which are strictly corresponding with the moisture change would provide the 
most accurate and reliable slope stability results and provide more reliable slope stabilization 
solutions. 
Keywords: Finite element (FE), Unsaturated soil, Soil Water Characteristic Curve (SWCC), Soil-
moisture diffusivity, Slope stability. 
1. Introduction 
The exploitation and processing of mineral 
resources plays an important role in Vietnam’s 
economy. Since the year of 2000, it has 
contributed 9.6-10.6% of the country's GDP 
(Tran Trung Kien and Pham Quang Tu, 2011) 
[1]. However, in reality, the industry has caused 
unexpectedly serious adverse consequences to 
the environment. Two projects on bauxite 
mining and alumina processing in the Central 
Highlands of Vietnam started in 2009, i.e., Nhan 
Co bauxite project in Dak Nong province and 
Tan Rai project in Lam Dong province. 
Environmental incidents in the two project sites 
have occurred continuously. 
An overflow of a chemical fluid from the 
workshop on chemical mixing workshop in 
August 2011 (dantri.com.vn, 2011) [2] and a 
leakage of a chemical fluid from a chemical fluid 
pipe into the surrounding area in February 2016 
(tuoitre.vn, 2016) [3] polluted some fish farms 
and the groundwater. On the 8th October 2014, 
an incident occurred with the basin of bauxite 
tailings in Tan Rai bauxite factory in which 
about 5,000 m3 of sludge was spilled out into Cai 
Bang lake, a 10-million-m3 irrigation reservoir 
which supplies water for the Tan Rai alumina 
factory. On the 23rd July 2016, an overflow of a 
chemical fluid occurred during the testing 
operation of Nhan Co alumina factory, as a 
consequence, 9.58 m3 of the chemical fluid was 
partly spilled into the Dak Yao stream (Figure 1). 
Figure 1. Nhan Co alumina factory, bauxite mine, red mud basins and landslide area. 
N.V. Hoang et al. / VNU Journal of Science: Earth and Environmental Sciences, Vol. 37, No. 1 (2021) 108-120 105 
Seven red mud basins (numbering from 1 to 
7) (Figure 1) have been designed for the Nhan 
Co alumina factory (Vinacomin, 2013) [4] with 
the total volume of 17.3 millions of cubic meters. 
At the present, only the red mud basin No. 1 with 
capacity of 0.26 millions of cubic meters is under 
operation. However, since August 2017 a series 
of landslides has been occurring in the 
neighboring farming land (Figure 2) along the 
red mud basin No. 2 during 2017 - 2019 rainy 
seasons, which lead to a complete destruction of 
transportation road along the red mud basin 
(Figure 3). In the near future when the red mud 
basin No. 2 and others are launched into 
operation, the landslides would definitely cause 
a possible overflow of the highly toxic fluid to 
the surrounding environment. 
Figure 2. 2018 landslide in farming slope in the red 
mud basin No. 2 - Nhan Co alumina factory. 
Figure 3. 2019 landslide on a road along the red mud 
basin No. 2 - Nhan Co alumina factory. 
As rain occurs, soil moisture increases 
because of the infiltration of the rainwater. The 
increased soil moisture leads to the decrease in 
the soil strength, i.e., cohesion and angle of 
internal  ...  and the bottom of the profile, 
while the suctions were determined through 
VWC and were determined by Eq. 15 with the 
given SWCC parameters. The initial VWC and 
suctions are presented in Figure 10. The VWC 
along the soil profile is presented in the upper 
abscissa versus the depth in the left ordinate. The 
suction is presented in the right ordinate versus 
VWC in the upper abscissa (the VWC and 
suction at the depth of 5.5 m are illustrated: the 
horizontal line going through the depth of 5 m 
and crossing the VWC curve gives VWC of 
0.508 on the upper abscissa, the vertical line 
going through VWC of 0.508 crossing the 
suction line gives the suction of 16.4 kPa on the 
right ordinate). 
Figure 9. Typical water content profiles for different 
soils [25]. 
Figure 10. Field initial water content and suction 
profiles for the modelled soil. 
4.3. Boundary conditions 
The first-type of boundary condition is 
assigned to the upper model domain of the 
ground surface since the moisture transfer model 
is to be carried out during excessive rainfall 
events. That is, the soil moisture content at the 
ground surface is corresponding to the saturated 
water content. As previously mentioned, during 
Aug. and Sep. in 2018 and 2019 a series of 
landslides occurred as a consequence of 
continuously heavy rain. During the 92 days 
from 1st Jul. to 30th Sep., there were 92 and 88 
rainy days in 2018 and 2019, respectively 
(Figure 11 and 12) (Dak Nong meteorological 
station) [26]. It means that it rained almost every 
day. The average daily rainfall for that period 
was 18.2 mm and 16.2 mm in 2018 and 2019, 
respectively. 
The lower model domain is also set up as the 
first-type boundary condition since it directly 
contacts with below-lying Upper Holocene 
aquifer. This is also in accordance to many other 
authors who assigned the boundary as either 
fixed water table, free drainage (unit gradient), 
or head dependent (refer to Carrera-Hernández 
et al., 2012) [27]. 
N.V. Hoang et al. / VNU Journal of Science: Earth and Environmental Sciences, Vol. 37, No. 1 (2021) 103-115 
112 
Figure 11. Jul-Sep daily rainfall in 2018 and 2019. 
Figure 12. Jul-Sep daily rainfall in 2018 and 2019. 
4.4. Simulation results of spatial and temporal 
soil moisture distribution 
The output of the moisture simulation is in 
the form of moisture content and suction along 
the depth for a given time interval, e.g., in hours 
or days. For the purpose of description, a soil 
moisture distribution curve at a concretely given 
time is used, or exactly at the 5th day from the 
beginning of the moisture transfer (Figure 13). 
For a concretely given time, the moisture 
distribution curve has three distinguished parts: 
i) an upper part, from ground surface to the depth 
Z1, with a saturated moisture, ii) an intermediate 
part, from the depth Z1 to the depth Z2, with the 
moisture in between the saturated moisture and a 
maximal field moisture at Z2 (C0,Z2), and iii) the 
field moisture from the depth Z2 and deeper. In 
the geotechnical analyses such as infiltration 
deformation and slope stability analyses, the 
intermediate part of the soil moisture distribution 
curve is of a special consideration due to the 
change of the shear strength parameters due to 
the change of the soil moisture. The upper part 
of the soil moisture curve has a saturated 
moisture and would have the shear strength 
parameters which are corresponding to the 
saturated soil. The lower part of the soil moisture 
curve has the natural field moisture, the shear 
strength parameters of which are 
correspondingly determined for during the 
geotechnical investigation 
The spatial and temporal soil moisture 
distribution curves for the whole simulation time 
of 90 days are presented in Figure 15. 
From Figures 13&14 it can be seen that the 
intermediate part length (from the depth Z1 to the 
depth Z2) increasingly changes with the time of 
the rainwater infiltration: at the 5th day the 
length is around 0.82 m and at the 30th day is 
around 1.2 m (Figure 14). Therefore, for a 
particular time since the beginning of the 
moisture transfer due to the rainwater 
infiltration, a certain length of the intermediate 
part needs to be specified with moisture for 
estimation of the soil shear strength parameters 
for the intended geotechnical analysis. 
Figure 13. VWC distribution curves for a concretely 
time at the 5th day. 
0
10
20
30
40
50
60
70
80
90
100
7/1 7/11 7/21 7/31 8/10 8/20 8/30 9/9 9/19 9/29
D
a
y
ly
 r
a
in
fa
ll
 (
m
m
)
Month/day (2019)
0.0
0.2
0.4
0.6
0.8
1.0
1.2
1.4
1.6
1.8
2.0
0.48 0.49 0.50 0.51 0.52 0.53 0.54 0.55 0.56
D
e
p
th
 (
m
)
Volumetric water content
Initial
After 5 days
0,Z2 s
Z1
Z2
N.V. Hoang et al. / VNU Journal of Science: Earth and Environmental Sciences, Vol. 37, No. 1 (2021) 103-115 
113 
Figure 14. Temporal VWC distribution curves 
during 90 days. 
Figure 15. Wetting front curves in the slope along 
red mud basin No. 2. 
Figure 16. Typical slip circles in regard to soil 
moisture and groundwater. 
Figure 16 presents three typical circles in the 
slope safety analysis: circle 1 is entirely within 
the saturated which has lowest shear strength 
parameters; circle 2 has a long part in the soil 
mass with a natural field moisture, and circle 
has a long part in the saturated soil mass with 
a buoyancy force effect. Among the three 
circles, circles 1 and 3 would have lowest 
factor of safety. 
5. Concluding remarks 
The finite element modeling of simulation of 
moisture transfer in unsaturated soils is a 
complicated and powerful tool which is capable 
of accurately determining the soil moisture 
distribution in time and space. The SWCC 
parameters are also time-consuming and 
expensively determined in the laboratory. Those 
two aspects in combination are challenging the 
unsaturated soil mechanics and hydraulics. 
Those laboratory test experiments and numerical 
modeling have been carried out for the slope of 
the red mud basin No. 2 in Nhan Co alumina 
factory in Dak Nong province. 
The SWCC parameters of the study soil are 
a=60, n=1.10 and m=1.05, which have been 
determined by the best fitting using the least 
squared error method, the mean squared error of 
which is very low and equal to 0.00011. The 
optimal finite element size for the soil under 
consideration is around 0.03 m, which had been 
selected based both on the analytical error 
estimation and on the numerical simulations 
with different element sizes. The results of the 
application to the study area's slope has shown 
that the wetting front depth can be up to 8 meters 
for 90 days of moisture transfer in the due to the 
rainfall. The wetting front depth and the length 
of the intermediate part of the moisture 
distribution curve have increased with the 
infiltration time. The soil moisture distribution 
with a depth is an essential information to have 
soil strength parameters for the slope stability 
analyses. Once, the moisture distribution at a 
particular time since the beginning of the 
0
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
0.48 0.49 0.50 0.51 0.52 0.53 0.54 0.55 0.56
D
e
p
th
 (
m
)
Volumetric water content
N.V. Hoang et al. / VNU Journal of Science: Earth and Environmental Sciences, Vol. 37, No. 1 (2021) 103-115 
114 
moisture transfer due to the rain, the 
intermediate part of depth needs to be specified 
with moisture content for the estimation of the 
soil shear strength parameters for the intended 
geotechnical slope analysis. The slope stability 
analysis with the soil shear strength parameters 
which are strictly corresponding with the 
moisture change would provide the most 
accurate and reliable slope stability results and 
provide more reliable slope stabilization 
solutions. 
There is a limitation matter remains the the 
study and would be interesting to be studied in 
derails. The upper boundary condition was set to 
be the first kind, i.e., a specified constant VWC 
during the entire modelled time domain based on 
the reality of 92 days of continuous rainfall, the 
mean rainfall intensity of which is 16.2 mm  
18.2 mm. The evaporation and excess moisture 
in the form of condensation for the air 
temperature drops below its dew point at the 
ground surface need to be identified for more 
precise specification of the boundary condition 
and boundary value of the ground surface. In this 
case, the upper boundary condition would be 
specified variable VWC or specified variable 
moisture flux etc. The effect of the boundary 
condition and boundary value on the moisture 
transfer results would assist in selecting a more 
appropriate boundary condition of the reality, 
Acknowledgement 
This paper was completed with the support 
of the project Code KC.08.23/16-20 financed by 
Ministry of Science and Technology of Vietnam. 
The authors would like to express a sincere 
gratitude for the support. 
References 
[1] T.T. Kien and P.Q. Tu. The role of NGOs in policy 
criticism: The case of the Consultancy on 
Development Institute in the sustainable 
management of mineral resources in Vietnam. 
Report of the Consultancy on Development 
Institute (CODE). Ha Noi, Vietnam (2011) (in 
Vietnamese). 
[2] https://dantri.com.vn/kinh-doanh/boxit-tan-rai-
tran-hoa-chat-ra-moi-truong-su-co-nho-thiet-hai-
lon-1317152926.htm (archived on 24 Sep. 2011) 
(in Vietnamese). 
[3] https://tuoitre.vn/nha-may-boxit-tan-rai-ong-chua-
xut-vo-do-lao-hoa-khop-noi-1052089.htm 
(archived on 16 Feb. 2016) (in Vietnamese). 
[4] Vinacomin. Environmental impact assessment for 
Nhan Co - Dak Nong bauxite mining. Dak Nong, 
Vietnam (2013) (in Vietnamese). 
[5] L. Rober, Schuster and R.J. Krizek. Landslides, 
analysis and control. National Academy of Sciences, 
Washington (1978); pp. 234. 
[6] K. Terzaghi, R. B. Peck and G. Mesri. Soil 
mechanics in engineering practice. John Wiley & 
Sons Inc. (1996). pp. 512. 
[7] Z. Huzhu, L.Hanbing, W.Jing and D. Weizhi. 
Investigation of the effect of water content and 
degree of compaction on the shear strength of clay 
soil material. Functional Materials. 24(2) (2017), 
290-297. ISSN 1027-5495. 
[8] Pelageia Iakovlevna Polubarinova-Kochina. 
Theory of ground water movement. Publishers 
"Science" (in Russian) (1977) pp. 664. 
[9] L.A. Richards Capillary conduction of liquids through 
porous mediums, Physics, 1 (1931) 318-333. 
[10] J.R. Philip Theory of infiltration. In Advances in 
Hydroscience Ven Te Chow (editor). Volume 5, 
Academic Press, New York (1969) 215-296. 
[11] B.H. Gilding Qualitative Mathematical Analysis of 
the Richards Equation. Transport in Porous Media 
5(1991) 651-666. 
[12] O.C. Zienkiewics and K. Morgan Finite Elements and 
Approximation. John Willey & Sons (1983) pp. 328. 
[13] W.R. Gardner. Some steady state solutions of the 
unsaturated moisture flow equation with 
application to evaporation from a water table. Soil 
Science, 85(1958) 228-232. 
doi:10.1097/00010694-195804000-00006. 
[14] D.G. Fredlund and A. Xing. Equations for the soil-
water characteristic curve. Canadian Geotechnical 
Journal, 31(3) (1994) 521-532. 
[15] D.G. Fredlund., A. Xing and S. Huang Predicting 
the permeability function for unsaturated soils 
using the soil–water characteristic curve. Canadian 
Geotechnical Journal, 31(3) (1994) 533–546. 
[16] E.C. Leong and R.H. Permeability functions for 
unsaturated soils. Journal of geotechnical and 
geoenvironmental engineering. Dec. 1997. (1997) 
1118-1126. 
N.V. Hoang et al. / VNU Journal of Science: Earth and Environmental Sciences, Vol. 37, No. 1 (2021) 103-115 
115 
[17] N.V. Hoang. Study on building finite element 
modeling software to simulate groundwater flow 
and pollutant and salt transport by groundwater-
application to the Vietnam Central coast plain. 
Vietnam NAFOSTED and Vietnam Ministry of 
Science and Technology. Project's code: 
ĐT.NCCB-ĐHƯD.2012-G/04 (2018) (in 
Vietnamese). 
[18] Nhan Co alumina joint stock company. 
Environmental impact assessment for Nhan Co 
alumina factory in Nhan Co commune, Dak R'lap 
district, Dak Nong province. Dak Nong, Vietnam 
(2009) (in Vietnamese). 
[19] Thuyloi University, Hanoi.  
[20] TCVN 8723 : 2012. Soil for hydraulic engineering 
construction - Laboratory test method for 
determination of permeability coefficient of soil 
(2012) (in Vietnamese). 
[21] ASTM D6836-02: Standard Test Method for 
Determination of the Soil Water Characteristic 
Curve for Desorption Using Hanging Column, 
Pressure Extractor, Chilled Mirror Hygrometer, 
and/or Centrifuge. Annual Book of ASTM 
Standards, Volume 04.09 (2002). 
[22] N.T.N. Huong and T.M. Thu. Determining shear 
strength of an unsaturated soil by the direct shear 
tests. J. of Water Resources & Environmental 
Engineering 42(9) (2013) (in Vietnamese). 
[23] C.W. Fetter. Applied Hydrogeology. Prentice Hall-
Upper Saddle River, NJ 07458. (2001) pp. 598. 
[24] P.S. Huyakorn and G. F Pinder. Computational 
method in subsurface flow. Academic Press Inc. 
(1983) pp. 473. 
[25] J Bear. Computer-Mediated Distance Learning 
Course on modeling groundwater flow and 
contaminant transport. Topic D: modeling flow in 
the unsaturated zone (2000). Faculty of Civil 
Engineering Technion-Israel Institute of 
Technology. Haifa 32000, Israel. 
[26] Dak Nong meteorology station. Monitored daily 
precipitation of the 2018 and 2019 years. Dak 
Nong, Vietnam (2020) (in Vietnamese). 
[27] J.J. Carrera-Hernández, B.D. Smerdon and C.A. 
Mendoza. Estimating groundwater recharge 
through unsaturated flow modeling: Sensitivity to 
boundary conditions and vertical discretization. 
Journal of Hydrology 452–453 (2012) 90–101.