Prediction of Groundwater Level Fluctuation towards Rainfall Induced Landslide: Case of Blue Nile Gorge, Central Ethiopia ()

The groundwater level fluctuation was studied in a complex geological
setting region where a frequent landslide is observed in a rainy season. Steady
and transient states of flow are modeled with different hydrogeological parameters.
The models are calibrated to satisfy the observed field conditions and expected
results from the scientific point of view. The results reveal that the groundwater level
fluctuation and flow direction in the region are complex. In limited areas, the
fluctuation of groundwater is significant from season to season while in others
the level remains stable in all seasons of the year. Following that, the result
of groundwater flow model was exported to GeoStudio to simulate the slope
stability of selected slope. The factor of safety was calculated using Slope/W.
The effect of pore-water on the factor
of safety was cross-checked by remodeling the slope without water. The results
and sensitivity analysis of slope stability confirm that the rise of
groundwater level decreases the factor of safety significantly only on critical
slope section.

KEYWORDS

Cite this paper

Yifru, B. and Ayehu, F. (2017) Prediction of Groundwater Level Fluctuation towards Rainfall Induced Landslide: Case of Blue Nile Gorge, Central Ethiopia. *Open Journal of Modern Hydrology*, **7**, 274-297. doi: 10.4236/ojmh.2017.74016.

1. Introduction

Regardless of improvements of science and technology in prediction and mitigating measures, landslides are still exact an economic and environmental toll in mountainous regions of Ethiopia. These hazards are becoming serious concerns to the public and to the planners and decision makers at various levels of the government. However, so far, little efforts have been made to reduce losses from such hazards [1] . The reason behind this is partly complexity of the processes driving slope failures and our inadequate familiarity of the underlying mechanisms. Numbers of factors have been investigated in the understanding of the driving force of landslides in the country much more in the Blue Nile Gorge (BNG). According to [1] , even though, many factors contribute to slope failure most of the stability problems in the country are from rainfall infiltration: earthquake triggered landslides are negligible.

BNG, along Gohation-Dejen road, is a well-known area in the country for frequent slope instability incidents. In the area, it is very common to see slope failure events that delay transportation during or right after a rainy period [2] [3] . Massive columnar jointed basalt, groundwater, uncontrolled surface runoff, joints of rocks, and the abundance of marl and shale within hard rocks are the leading reasons for slope failure. More recently, this slope instability problem has been damaging the road sections, bridges, and tillage [2] . In the region, still, there are active landslides [2] .

There are ongoing and completed studies in the area to investigate the problem. The recent studies by Japan International Cooperation Agency (JICA) and Geological Survey of Ethiopia (GSE) identified landslide areas in the region are in continuous movement. Some of these are up to two kilometers wide, putting into jeopardy this vital link [4] . Most of the land features in which the road passes through are mountainous, unstable, and highly susceptible to sliding. Even though several kinds of literature conclude, major triggering factor for slope instability in this landscape is groundwater level increase from infiltration into the soil, there is no comprehensive groundwater study in the area, except the above-mentioned field investigation for the purpose of the road protection. Nearly all researchers tried to see the failure from the landscape and geological perspective. Nevertheless, extreme hydrological event (groundwater flow) is a determinant factor in such landscape on deferent phenomena [5] .

Groundwater flow is dynamically connected to the hydrological cycle through different recharge processes. As part of the hydrological cycle, it is continuously in motion from regions of recharge to discharge sites. Landslip areas are the site of local groundwater storage and moisture retention in the region that would or else lose water quickly to neighboring valley owing to high slope. Depending on materials, landslide bodies hold a significant quantity of groundwater, which discharges mostly as diffuse flow or as focused discharge to springs. A typicalexample of affiliation of landslides with springs is that taking place within the BNG [6] .

Moreover, it is widely recognized that water is one of the foremost initiates of landslip commonly in the highlands. Considerable landslide studies had conducted so far everywhere in the world. Utmost of these researchers discuss results that water may have on slope stability such as lowering suction, rising groundwater stage and causing an upsurge in pore-water stress, seepage erosion, hydraulic uplift strain from underneath the landslide, and influence of water at the plasticity of the landslide. Furthermore, studies factor out that long lasting rainfall has an enormous effect on slope balance. A normal lower in precipitation results in a decreasing of the water level, as well as a diminution in the weight of the soil mass. Then again, an increase in precipitation will improve the extent of the groundwater, lesser shear trength, growth the load of the soil mass, and might increase erosion.

More recently, the development of computer technology and packages make models powerful equipment for environmental safety: sustaining the groundwater equilibrium, and protecting against intense land subsidence [7] . A model is a representation of a real system or process [8] . The result of a model may be inferred expectantly outside its database in comparison to empirical techniques wherein the problem is well described.

There are different forms of models, inclusive of the conceptual model, numerical or analytical model, and physical model. From these, conceptual and numerical models are widely in use in field problem analysis. The first one is a hypothesis for the way a system or process function. Whereas mathematical model is an illustration of a conceptual model with objects, forces, and events replaced with mathematical language.

Mathematical models usually solve partial differential equations representing conservation of mass, momentum, and energy. Models use diverse numerical procedures to solve the governing equations conditional on suitable initial and boundary conditions. Applications of a numerical set of rules that solve the governing equations are in a computer code that may be considered a generic model.

Even though numerical techniques are tools, which can resource in studying groundwater problem and might assist to increase our expertise of groundwater systems [9] , they need complicated computation and require more data of boundary conditions and material properties. Nonetheless, numerical methods can give more accurate results in comparison to analytical methods [10] .

The purpose of this study is to characterize groundwater stage fluctuations, in particular, previously instrumented landslides and adjoining slopes using coupled hydrologic-slope stability models.

2. Materials and Methods

2.1. Features of the Study Area

The study area is located in Blue Nile Basin (BNB), which is situated in the Northwestern part of Ethiopian. It is around 189 to 229 km from the capital city of the country, Addis Ababa. Highly important Addis Ababa-Debremarkos-Bahir Dar-Gondar-Metema-Sudan and Gondar-Tigray road that links North central and Northwestern part of the country with the capital and port of Sudan crosses this region (Figure 1).

The BNG landform is that of a basaltic lava plateau (Eocene period flood lava), with an elevation of about 2500 m, underlain by Mesozoic sedimentary

(a) (b)

Figure 1. The setting of study area: major river basins of Ethiopia (a), Gohatsion-Dejen area (b).

rocks of various origin and type as low as 1028 m above mean sea level [4] . Lateral slopes of BNG consist of several levels including cliffs and colluvial slopes. These gentle slopes develop on the areas of limestone and shale, which are covered by residual soil and colluvial deposits. Although the sedimentary and volcanic rocks in the area are exposed largely as symmetrical stratigraphy on both sides of the Blue Nile River (BNR), the detailed sequences are unevenly distributed. The sequence in the area is not disturbed due to major faults and is generally horizontally stratified. Figure 2 shows a schematic geological cross segment of the BNG [11] .

The climate in BNG is warm and temperate. The average temperature is 19.5˚C with a maximum of 37.7˚C. The rainy season is from the end of May to September, with July and August accounting for approximately half of annual precipitation. The average, minimum, and maximum rainfall in the area are

Figure 2. Schematic geological sections of the BNG (Ayalew L and Yamagishi H, 2003) with few modifications.

1270, 1209 and 1362 mm per year respectively over 22 years record from 1994 to 2015. Figure 3 shows twenty-two years average monthly rainfall pattern of the area. The precipitation in the region is fairly regular year to year. At the site of interest, there are five rainfall-measuring stations namely Abay Sheleko, Church, Bridge, Filiklik, and Gohatsion. For groundwater recharge from rainfall adjustment, the coverage of the area by each station is approximated by Thiessen polygon method as shown in Table 1 and Figure 4 below. This result is going to be used later in the adjustment of recharge in groundwater flow analysis. Groundwater recharge in this basin is 100 mm to 303 mm per year [12] . The potential evapotranspiration of the area is also shown in Figure 5.

The BNR crosses this area. It is the largest and perennial river in this site of interest. There are also small streams, which run from the hilly side and dry out after a few period of the rainfall season. The runoff and recharge in the basin, in general, depend on the topography and rainfall amount. Relatively flat area of the region is dominantly covered with agricultural crops while the hilly is covered with lightly dispersed small trees.

2.2. Data Processing and Model Setup

The objective of the research achieved via a serious of tasks systematically. The necessary data are collected both from organizations and from the site. The collected data include daily rainfall record, hourly ground water head, hydrogeological characters of the aquifer, and geo-mechanical properties. All the data were collected from Ethiopian Road Authority (ERA), Metrological Agency of Ethiopia (MAOE), and Geological Survey of Ethiopia (GSE) except soil strength properties. These data were examined, organized, and prepared as the form that is acceptable to the software using Microsoft Excel sheet, Surfer, and Global Mapper applications. Following that, all data organizations and modeling are done accordingly. Brief systematics of methods followed is depicted in Figure 6.

Figure 3. Monthly average rainfall (mm) of the study area over a period from 1994 to 2015.

Table 1. Average areal daily rainfall of the study area.

2.2.1. Groundwater Flow Modeling

The selection of the ideal method for modeling a certain problem relies upon on many aspects, which includes the required degree of precision, modeling expenses, the simplicity of use and understanding of the results [13] . MODFLOW, finite-difference model, was used to simulate groundwater flow. It is broadly used modular 3-dimensional block-centered code [13] in Processing MODFLOW modeling environment [14] . Processing MODFLOW is a graphical interface that integrates MODFLOW with numerous programs to simulate a selected characteristic of a hydrologic system. It enables representing an extensive range of drainage conditions, geometry, configurations, and different hydraulic settings. It gives a variety of boundary conditions, such as constant head, recharge, evapotranspiration, drains, rivers, etc. [15] .

An aquifer might be confined or unconfined, or set as convertible from confined to unconfined and the other way round. In this extent, non-uniform and the highly variable vertical relationship between the aquifer layers and the overlying confining beds makes the assumption of different model layers difficult. Therefore, the region was simulated in steady state and unsteady state,

(a)(b)

Figure 4. Thiessen polygons (a) and Monthly Rainfall at the five selected stations (b).

Figure 5. Potential Evapotranspiration (PET) of BNG.

Figure 6. Systematics of the research.

unconfined aquifer, and single layer approach, which is more practical to avoid complications because of many unknown parameters and geometrical condition of multilayer approach.

The basic principle of groundwater flow fundamentally lies on Darcy’s law. When this law put together with continuity principle it can be described by the general flow equation in 3-dimensions for a heterogeneous anisotropic material. The partial-differential equation of groundwater flow used in MODFLOW is [16] :

$\frac{\partial}{\partial x}\left({K}_{xx}\frac{\partial h}{\partial x}\right)+\frac{\partial}{\partial y}\left({K}_{yy}\frac{\partial h}{\partial y}\right)+\frac{\partial}{\partial y}\left({K}_{zz}\frac{\partial h}{\partial z}\right)\pm W={S}_{s}\frac{\partial h}{\partial t}$ (1)

where: K_{xx}, k_{yy} and K_{zz} are values of hydraulic conductivity in the x, y and z directions along Cartesian coordinates Axes, h is hydraulic head, W is volumetric flux per unit volume and represents sink and/or sources, S_{s} is the specific storage of the porous material and t is time.

MODFLOW solves Equation (1) using the finite-difference approaches where in the groundwater flow system is divided right into a grid of cells.

1) Spatial Discretization

Model grids discrete the continual natural system into cells that enable to do the numerical simulation. The distance from node to node, that is named grid resolution, ought to be receptive to sharp changes. The overall size of the grid ought to be capable define the problem and therefore the results of the procedure consistent with modeling purpose, however not so large to cause excessive preparation and computation needs.

In this study, the region is replaced with a set of discrete nodes in a grid pattern covering the modeled area. The grid includes 324 rows and 292 columns overlaying on the 86.431 square kilometers. Moreover, it is refined around the observation wells to be able to increase the calculation precision of the model. Since the aquifer shape is not square, it shows that the model domain in MODFLOW exceeds the study region defined in the conceptual model. Consequently, we defined the IBOUND by first setting up the lateral extent of the formation using the catchment boundary map and assigned a cell as active if the formation enclosed greater than 50 percent of the cell area.

2) Boundary Conditions, Model Parameters, and Stresses

Following model area discretization the compulsory information for each cell was imported, such as top and bottom layer of the aquifer, hydraulic conductivity, specific yield, and porosity of the formation and other trivial. The data for the model maybe classified as spatial and historical. The spatial input comprises aquifer features, such as boundaries, hydraulic conductivity, location of wells, recharge area, drainage area, etc., whereas the historical input includes time dependent records especially for transient simulation.

The boundary condition, grid dimensions, initial aquifer properties, and time step features are specified as the basic need of the model. The surface topography was derived from a survey carried out over the area and the aquifer extent was structured based on the available details of the boreholes stratigraphic entity. The mean thickness of unconfined aquifer layer was 400 m comprising different hydro-geological fields.

The top of the aquifer layer is assumed on average 12 m below the ground elevation and drawn from the 30 m by 30 m resolution digital elevation model (DEM). The DEM data is processed to generate grid file to be in a format compatible with the Processing MODFLOW Pro (Version 8.0.15) and imported to the model with elevation referenced according to its geographic position. The imported grid elevation has an elevation range of 1028 to 2500 m. In order to prevent drying of cells during modeling, elevated zones were given different thickness based on their depth from the correctional view, relatively higher thickness at higher cells (Figure 7).

Hydraulic boundary conditions are presentation of dependent variable (head) or the derivative of the dependent variable (flux) at the boundaries of the problem domain. Setting the boundary is the critical step on the modeling. Three types of boundary conditions were used to define the groundwater flow system in the study area: no-flow boundaries, specified-flux boundaries, and constant head boundaries. Geologic or hydrologic barriers to groundwater flow were simulated using no-flow boundaries (Figure 8).

The major perennial river (BNR) is another significant boundary condition, which facilitates the movement of groundwater to and from the river. The flux of water between groundwater systems and rivers is generally dependent on the hydraulic head in the aquifer and is simulated as a head-dependent flux boundary. In order to simulate the interaction between river and groundwater, the model used the river package [16] . Estimated riverbed conductance was based on model calibration.

The effect of springs was modeled using drain package, which removes groundwater from the aquifer at a rate proportionate to the head difference between the aquifer and the springs i.e. when the hydraulic head in the aquifer is greater than the drain (spring) elevation, groundwater flows towards the drain.

The very indispensable parameter in the aquifer system is the hydraulic conductivity that defines the flow rate of the groundwater in the aquifer system. It

(a)(b)

Figure 7. Top (a) and bottom (b) layers of the aquifer.

can be defined as a unit of discharge of water through a unit area of a porous medium under unit hydraulic gradient measured at right angles to the direction of flow. Hydraulic conductivity is a function of both the medium and the fluid. The model uses the spatial distribution of the hydraulic map described in the conceptual model to begin the model simulation. In this model, groundwater flow within the layer was assumed horizontal (Figure 9).

Recharge is a specific flux boundary, which is independent of the head of the cell, but MODFLOW consider it as a property for spatially distributed all over the model area. Recharge to the model consisted of infiltration from direct precipitation, stream infiltration draining the Northeast and Southeast escarpments.

Figure 8. Location of observation wells and springs, model boundaries.

Figure 9. Geological zones used for hydraulic conductivity estimation.

Recharge was applied to the active model area as a spatially varying, specified flux to the highest active cell. In general, precipitation recharge varies spatially with land surface permeability, which is a function of soil characteristics and land use, and spatial distribution and intensity of rainfall. For the year period (365 days), different recharges for the five recharge zones as input into the model. The recharge is applied with recharge package. As it is a boundary property, the recharge has a mathematical code that assumes the volumetric rate of flow into cell described with multiplication of the recharge rate by the horizontal area of the cell (Figure 10).

2.2.2. Steady State Model Development and Simulation

Steady state flow is the state in which flow parameters passing a given point per unit of time remains constant. An initial steady state condition is required for transient state modeling of groundwater flow. The required data including recharge, initial hydraulic head, evapotranspiration, hydraulic conductivity etc. were imported to the model system. After importing the necessary information and running the model, the hydraulic heads are the main results. To simulate steady state groundwater flow initial heads are required as starting values. Measured groundwater levels, collected from the observation wells in the study area, were used to retrieve the initial hydraulic head of the aquifer. These values were used to the model as initial groundwater heads.

Eleven observation Piezo-meters and selected springs, where the water level measurements were averaged to daily, were used for model simulation. The wells have different head value within a short distance difference between them. In addition to that, some wells show comparatively significant difference in a year whilst others have approximately same groundwater level in all months of the year. Nevertheless, all measuring wells show a continuous increase in level from 2012 to 2014 except wells labeled as B2823 and B0513. Model simulations were completed over a daily period (Figure 11 and Table 2).

2.2.3. Transient Model Development and Simulation

Transient groundwater flow is a dynamic system, in which flow parameters such as inflows, outflows, and groundwater storage vary with time. In transient

Figure 10. Average historical groundwater level below the earth surface (data from ERA).

Figure 11. Recharge zones and values in m/day.

Table 2. Water budget of the model at steady state flow simulation.

groundwater flow simulation, a steady initial condition, which is a result of steady state model, has been adopted. A specific yield of 0.468, a maximum evapotranspiration of 1.20E-08 m/d, and specific storage of 0.014 were used. One-year long simulation is divided into 12 stress periods; each stress period represents a month, while the interval of each stress period is divided into days. The total time steps, therefore, equals 12 months, while the total simulation time equals 365 days (Table 3).

In addition to this, a model maybe used to predict the upcoming hydraulic response of an aquifer, under different scenarios, such as decrease/increase in

Table 3. Average monthly groundwater recharge from rainfall.

recharge due to land use change. It was tried to forecast groundwater level in the coming five years by dividing the stress period into 60 stress periods i.e. one-year simulation is divided into 12-period lengths. Variable groundwater levels were noticed at different locations of the aquifer but the variation of the level with time is very slow and minor (for example compare Figure 12 and Figure 13) even ignorable (Table 4).

2.3. Slope Stability Modeling

This part of the study focuses on limited landslide area where attempts were made to investigate the effect of groundwater level on the stability of a slope. The area is located on the boundary of the basalts. Highly weathered basalt is intercalated between the basalt and the limestone. It is sand or mud of soft particles prone to liquefying in the case of containing water. The layer is very soft and erosive, which would trigger slope collapses [17] . The slope is in very high landslide zone [18] .

The limestone that is fine to medium grained is widely distributed in the region. In addition to that, it is rich with Colluvial deposit, mainly composed of basalt boulders, gravel, sand, mud, and clay including black and organic soil as “black cotton soil”, which have come from the mountainous side cliff. Water permeability of the layer is relatively high. The deposit is a source of the rock falls and slope failures in the rainy season [4] (Figure 14).

Coupled modeling approaches increasingly become popular in slope analysis. This includes combined limit equilibrium stability analysis and groundwater flow (Hydrologic-Slope stability) modeling. In this paper, the analysis was performed for a single slip surface as a way failure back analysis using a module of GeoStudio 2012: Slope/W to calculate the factor of safety by Limit Equilibrium Method (LEM). Slope/W allows carrying out limit equilibrium slope stability analysis of natural and engineered slopes. The program has several methods such as Bishop’s Modified method, Janbu’s Simplified method, Spencer method, Morgenstern-Price method and others. These techniques are functional to circular, composite, and non-circular surfaces. LEM approaches are the most widely used for analyzing slope stability. The simplicity and versatility of the LEM rest with the concept that the geometry of the possible failure surface in a

Figure 12. Contour of computed groundwater level at steady state.

Figure 13. Contour of groundwater level at transient state for month of August (Stress period 8, time step 15).

Table 4. Water budget of the model at transient state of flow.

slope is identified and the slope can be discretized into ﬁnite vertical slices. Each slice is then analyzed using principles of force and/or moment equilibrium [19] for its role to the rigidity of a slope.

Figure 14. Landslide hazard zone map of BNG [18] .

In this work, Morgenstern-Price Method selected, which consider normal and tangential equilibrium including moment equilibrium for each slice in circular and non-circular slip surfaces [20] . It solves for the factor of safety using the summation of forces tangential and normal to the base of a slice and the summation of moments about the midpoint of the base of each slice. The equations were written for a slice of small breadth. The force and moment equilibrium equations were combined and a modified Newton-Raphson numerical technique [21] .

Moreover, Slope/W is the foremost software product for computing the factor of safety of earth and rock slopes [21] . It enables to analyze both simple and complex problems for a diverse of slip surface shapes, boundary conditions, soil properties, analysis techniques, and loading environments.

After the sensitivity analysis to assess the effects of mesh reﬁnement, apt mesh conﬁgurations were fixed aimed at this site as shown in Figure 15. It is an unstructured triangular finite-element mesh with the intention of practically all elements would be equilateral triangles, which are expected exceedingly improving computational precision.

The important parameters in the hydro-mechanical framework are the drained cohesion c’ and friction angle ϕ’ deﬁned in the Mohr-Coulomb failure criterion. Direct shear tests were conducted under saturated conditions to get shear strength parameters of soil. The parameters are also estimated by adopting inverse analysis [22] . The values for these hydro-mechanical properties are listed in the following Table 5.

As it is depicted in Figure 16 the relationship of the three important parameters under consideration seems are not in good agreement. However, if we

Figure 15. Gridded Slope geometry, material used and position of pore-water pressure line in the specific slope.

Figure 16. Time (from 2012 to 2014) trend of rainfall, groundwater level, and landslide at selected very high landslide zone of the region.

Table 5. Hydro-Mechanical Parameters.

examine the data over projected longer time (even in the period given) it is clear that the difference between extreme events is time. This is outlined in the introduction section as landslides are occurring at the end of rainfall seasons. In addition to that, control of topography is high in the area since the groundwater is shallow and direct rainfall infiltration into the aquifer is slower.

Most commonly, slope stability analysis result is measured with a factor of safety (FS) which is the ratio of resistance force against landslide soil mass to force when a landslide triggered along the slip surface. It enables to tell how close or far a slope is from failure. For example, FS is equal to one means resistance forces and sliding forces are balanced. If FS is, greater than one means the landslide is stable whilst FS less than one is unstable or sliding. A method of slices was used to compute the FS along failure surfaces and to search a critical slip surface (a surface with the lowest FS).

Initially, the effects of topography and geology were inspected using two types of materials without pore-water pressures. In this set-up, the least stable area is located on the steeper slope. The minimum FS obtained from the analysis is 1.226 (Figure 17).

For the second analysis pore-water, pressure field is taken from the result of MODFLOW and incorporated into the model. Realistic model is obtained by combining pore pressures with heterogeneous strength properties. The result shows the least stable area where pore pressure is locally elevated in and the slope is relatively steeper. Figure 18 show strength and water pressure plummet and upswing considerably at the different longitudinal profile of the slope geometry. The effect of pore-water development is analyzed and after many trials; it is noticed that it has the capability to reduce the factor of safety up to 18 to 22 percent in the particular site studied in general (considering the factor of safety of each slice). At critical slice, however, every less than half a meter increase

(a) (b)

Figure 17. Slope stability analysis results with no pore water pressure (a) and with pore water pressure imported from MODFLOW at steady state (b).

(a) (b)

Figure 18. Pore-water Pressure (a) and strength along the slope (b).

in groundwater level reduces the factor of safety nearly by 0.84 percent only. Results are compared with recorded landslides and there is acceptable correspondence.

3. Model Calibration and Sensitivity Analysis

Calibration is a procedure of modifying a model for a specific problem by using manipulating the input parameter, and initial or boundary conditions, within reasonable levels until the simulated model results intently fit the measured data. The typical approach is to select objective function that may be a measure of the agreement between observed and simulated facts, and this is directly or indirectly associated with the amendable parameters. Optimum parameters are obtained through minimizing the objective function. Minimization of the merit function (model calibration) may be done through trial and error or, as is becoming more famous, by the use of an automatic parameter approximation system. Calibrated model can reproduce facts within some subjectively satisfactory degree of accuracy.

In this study area, local information about the values of hydraulic conductivity is limited. It should be estimated based on the observed head and other hydrogeological parameters. Therefore, the area is zoned into twelve main lateral geological formations of the aquifer and the hydraulic conductivity values of these zones were defined as unknown parameters. On all zones, 1 × 10 − 7, 100, and 0.01 m/day were specified as the minimum, maximum, and initial values respectively in order to get optimized parameters using Parameter Estimation (PEST) model.

Regularized inverse problems are commonly solved using Parameter EST imation code PEST [23] to determine the value usually hydraulic conductivity or recharge by using historical data such as piezo-metric head or stream discharge. This model tries to estimate parameters by minimizing the residual. The objective functions are calculated as the sum of squared weighted residuals or as mean square error (MSE) [24] .

$MSE=\frac{1}{n}{\displaystyle {\sum}_{i=1}^{n}{\left({H}_{c}-{H}_{o}\right)}^{2}}$ (2)

where; n = number of data; H_{c} = computed value; and H_{o} = observed value.

In each iteration, to allow the user to decide sensitive and insensitive parameters, the PEST model also calculated parameter sensitivities. Many adjustments in the size of zones were made during the simulations in order to balance the sensitivities. Finally, a realistic fit was obtained between observed and computed piezo-metric heads. In Figure 19, measured head versus observed head is shown after calibration in the area in general. It shows good agreement between observed and computed values.

The model was calibrated in two steps, first using steady state simulations to estimate hydraulic conductivities of the model, then using a transient simulation to estimate the parameters that primarily affect fluctuations in flow. The model used a hydraulic conductivity within a range of 0.002 m/day to 19.82 m/day.

Figure 19. Measured head versus observed head.

(a) (b)

Figure 20. Calculated head versus observed head at transient state for stress period one and eight.

Adjusted values from the steady state modeling were then used in efforts to calibrate the model using transient settings, which included modifying stream conductance, recharge, evapotranspiration, and storage coefficients. Calibration of the transient simulation yielded estimates of storage coefficients and estimates for stream conductance. One-year imulation, without sinks, was divided into 12 stress periods; each stress period represents a month, wherein the length of each stress period is divided into days. The total time steps, therefore, equal 12 months, whereas the total simulation time equals 365 days. There were some verification data of two months. The existing data were used to calibrate the transient state of the model. The first step of the calibration was to assign an initial value for storage coefficients and specific yield for the model. Transient model calibration was done using trial and error procedure by changing the specific yield, storage coefficient, and with a very limited range of the hydraulic conductivity values. Acceptable performance of the model was observed through the transient simulation to fit between simulated and observed head in the observation wells (Figure 20).

4. Conclusions

The role of groundwater flow model is to define the equilibrium of inflow-outflow events so that changes in local groundwater flow rates and fluctuations in water level can be predicted. Regardless, in mountainous regions of Ethiopia, the scholarly work about groundwater is relatively low due to several reasons such as inadequate data, the complexity of aquifer geometry and inaccessibility of the area. In this paper groundwater table variation in the unconfined aquifer principally owing to rainfall infiltration, recharge from the surrounding, and interaction with the river was modeled.

The model is calibrated using historical data and we obtained are asonable variation of the groundwater table in response to the above-mentioned stresses. The results show that groundwater table is varying depending mainly on recharge from the hilly side of the extent. There is a contribution from rainfall infiltration too, but it is trivial, especially on Dejen side of the gorge. However, in a limited section of the region, the upsurge and drop of groundwater level during rainfall and dry seasons are noteworthy and the amount of water level difference in the area between the driest and wettest time is not significant except in limited regions. In addition to that, this work exhibited; though, in limited areas the recharge is high into the aquifer and may contribute to the slope instability at slopes that are classified as high landslide hazard zones, in the majority of the regional groundwater stays stable throughout the year, which implies that the rainfall infiltration in these sites has less effect on stability. Their geographical location, sensitivity analysis, and geological zone imply that they are receiving water from the northeast and southwest neighbor of the proposed model. Model results and field measurements show the geometry is in continuous movement. A slight rise in water level decreases the factor of safety of the slope substantially. Less than half a meter increase in water level, decrease the factor of safety up to approximately 0.84 percent. Furthermore, it is understandable that the pore-water development corresponds to the slope of the area.

From this study, we can conclude that in order to make the investigation of slope instability and mitigation measures coupling groundwater flow modeling and slope stabi