The Effect of Geometrical Features of Release Surfaces on the Stability of Tectonically Disturbed Deep Rock Slopes in an Albite Open Pit Mine

Abstract

An albite deposit formed in orthogneiss has been operated using open pit mining process since 1996. However, landslides have occurred following the excavation works and heavy rainfalls. In the current state, the slopes are still susceptible to mass movement. For this reason, this study aimed to determine the root causes that are influencing continuous movements. In this context, two geometric features of the strike-slip faults that form tectonic slices in the pit were proposed. The first one (α-angle) is the acute angle between the dip direction of the overall slope and the general trend of the strike-slip fault and it was used for the kinematic evaluation of the sliding potential of the tectonic slices. The validity of α-angle was confirmed by investigating its relationship with the factor of safety values. The second one (undulation angle, i) was used to determine the effective surface friction angle (ϕ_eff). Thereafter, ϕ_eff value was used in the deterministic finite element (FE) slope stability analyses as well as back analyses. Its effect on the stability of overall slopes was investigated for static and pseudo-static conditions with different water table levels. The FE results confirmed the monitored field displacement measurements. It was found that as α-angle decreases, the ability to resist sliding also decreases since the low α-angle promotes block formation. This study also demonstrated that using the undulation angles on mapped profiles of large-scale discontinuities (≥270 m) would be a promising approach for slope stability assessments of open pit excavations in tectonically disturbed metamorphic rocks.

1. Introduction

Geological structures such as folds, faults, and discontinuities are intrinsic parameters that play a critical role in the stability and behavior of both natural and engineered rock slopes [1,2,3,4,5,6]. When a rock mass shears along a structural plane, the shear strength of this plane is affected by the undulations and by the friction between the contact regions [7,8,9]. A considerable amount of literature has been published on the scale dependency of surface roughness and its effect on shear strength [8,10,11,12,13,14,15,16,17,18]. On the other hand, several researchers studied the effect of geometry such as multiple-faced slopes and concave and convex geometry of the slope facets on the stability by kinematic analyses and/or numerical modeling [18,19,20,21,22,23]. Hencher [24] proposed a method for quantifying the roughness of a joint as the deviation from the mean dip angle in direction of potential sliding at a scale of 40 cm. To date, only Antolini et al. [25] measured variations of the dip direction and the inclination of the bedding planes to consider large-scale (<70 m) undulation in numerical modeling. Recent developments in the field of monitoring facilitated preparing numerical models of 3D in-situ discontinuity surfaces with the aid of terrestrial laser scanners and image analyses to reveal the roughness characteristics [26,27,28]. However, a major problem with this kind of application is that undulation angles of the almost vertical surfaces cannot be directly measured due to high-dip angle and the application will be very challenging in the case of very long profiles. In addition, these studies have not dealt with geological structures in open pit scale (≥270 m) such as strike-slip faults that are dividing the overall slopes transversely and forming lateral boundaries for mass movements. Accordingly, undulation angle of the large geostructural features have not been quantitively treated as a geotechnical parameter in numerical modeling.

On the other hand, several researchers have investigated the types of failures; inherent and environmental causes of these failures occurring on the eastern slopes of the studied open pit in previous years. Tanyaş and Ulusay [29] examined planar and wedge-type failures on the benches with the aid of back analyses using the limit equilibrium theory. Kıncal [30] studied planar sliding on benches and its relationships between tension cracks on the upper slope face and discontinuity sets by stereographic projection technique. Kadakci [31] and Kadakci Koca and Koca [32] studied overall slope stability by numerical models considering 45 m deepening for ore production. Karagöz and Koca [33] studied the continuous displacements on the southeast slopes of the open pit using a GPS instrument. They also linked the fluctuations in monitored movements to the amount of rainfall. They revealed that the progressive stage of failure will impact six benches. Hydrogeological conditions were identified. Accordingly, remedial measures have been recommended. More recently, Karagöz et al. [34] studied the landslide that occurred on the eastern slopes in December 2018. They determined the extent and mode of the landslide. However, in 2021, the movements after water drainage measures (water pipes) still proceed. In addition, further production of ore material is planned, and the safe deepening of the pit can only be operated with a proper understanding of the special (unique) geological features in the field. For this reason, this study seeks to fill the gap in research on the scale of undulation by measuring the curvatures of the surface traces of the undulated, large (pit-scale, ≥270 m) structures on the map. Accordingly, the measured curvatures would be treated as undulation, hereby an undulation angle (i) representing each shear zone has been obtained. The concept of adding a small or block scale i-angle to the laboratory-measured internal friction angle (ϕ) to obtain the field shear strength was recommended by Hencher [24]. Therefore, the undulation angles were integrated into the two-dimensional (2D) FE modeling to simulate the shear strength of discontinuities. Additionally, another parameter (α-angle) representing the directional relationships between the shear zones and the slopes was investigated whether it contributes to slope instabilities or not in a comparative manner. In this context, wet zone determination, detailed fracture map, and mapping of landslide areas were thoroughly performed to associate the slope instabilities with the pit-scale shear zones (tectonic slices) as well as hydrogeological characteristics.

2. Geological Setting

The albite open pit mine is located in Aydın Province in western Turkey. In the open pit and its nearby, the tourmaline-rich orthogneisses are exposed which is the most typical lithology of the core series of Menderes Massif (Figure 1). The Menderes Massif was developed by compressional deformation [35]. It has an onion-shaped structure that is composed of gneisses in the core and covered by a Paleozoic schist envelope [36]. Core series in the study area are composed of medium-grained, albite-rich gneisses. They show a pronounced trend in the NE-SW direction (Figure 1). The mineralogical composition of leucocratic veins is albite (44–47%), orthoclase (1–5%), quartz (41–55%), and rutile/sphene (1–7%). They have a medium-grained (<0.4 mm) granoblastic texture. Ore bearing zone (with the mineralogical composition of Na-feldspar) was developed along an old compression zone.

There are three foliation planes and five joint sets in the open pit as follows; Foliation planes; F₁: 310–347/32–36, F₂: 273/40, F₃: 243/40, Joint sets: J₁: 332/80–90 (Shear joints), J₂: 190/70–90 (Shear joints), J₃: 250/88–90 (Shear joints), J₄: 280–290/74–90 (Compression joints), and J₅: 221/84 or 41/84 (Shear joints). On the other hand, nine strike-slip faults intersect the overall slopes of the open pit in the directions of NE-SW and NW-SE (Figure 1). These strike-slip faults are dividing the pit area into eight shear zones. Each of these zones represents a tectonic slice that is restricted by the nine strike-slip faults from north to south. In the NE of the open pit, the orientations of the strike-slip faults that perpendicularly intersect the normal fault-I are trending between N20E/90 and N86E/90. However, the trends of strike-slip faults are changing in the SE sector of the open pit as E-W/90 and N84W/90 (Figure 1). In this sector of the pit, the hanging wall also remains within the perimeter of the open pit slopes, which makes this sector very fractured and susceptible to mass movements.

Some of the strike-slip faults are differentiated from others due to the high lateral offset (>20 m). These tectonic slices contain many small (secondary) shear joints that were traced in the field elongating nearly parallel to the main shear joints and each has a different shape (Figure 1). These slices have the potential to cause stability problems owing to their geometric features. The strike-slip faults with high lateral offset which make it easier to slide, form the release surfaces of each tectonic slice. The value of lateral offset on the strike-slip fault which is restricting shear zone-6 is about 41 m. Therefore, the shear zone-6 coincides with the area of an old landslide that occurred on 19 December 2018 (Figure 2).

The old landslide developed almost parallel to the slope face (Figure 2). The failure surface was initiated from tension cracks developed from the discontinuities with the orientation of 280–290/74–90 planes and followed the foliation planes with the orientation of 310–347/32–36. Therefore, foliation planes are forming basal planes although they do not daylight on the overall slope face but locally on benches (α_bench > α_foliation ≥ α_overall). The strike-slip faults that intersect the foliation planes and bound the tectonic slices from the N-S borders worked as lateral release surfaces. On the other hand, almost vertical joint sets trending perpendicular to the strike-slip faults are the rear release surfaces of the landslides. Overall, the strike-slip faults are forming the lateral limits of slope failures for each tectonic slice. For this reason, the presence and degree of undulation angles of the surface traces of these faults will significantly affect the stability of slopes.

Hydrogeological Conditions

The orthogneiss rock mass is recharged by the leakage from the creeks behind and on the pit slopes. A preferential flow oriented parallel to the 32° to 40° inclined foliation planes is sometimes interrupted by the valley-parallel, steeply inclined joint sets (N-S oriented joint sets) and normal faults which act like water barriers. Accordingly, perched water levels can be observed in the winter. As the foliation planes are oriented towards the base of the pit and occasionally daylights on the bench slope faces, water seepage along the foliation planes is observed at several benches. The strike-slip faults also provide good drainage conduits for the eastern slopes along their strike. The 82°–84° inclined strike-slip faults that are vertically cutting the foliation planes restrict water flow from one shear zone to another and the water directly flows toward the pit base.

3. Geometrical Feature of the Strike-Slip Faults: α-Angle

The first measured geometric feature is α-angle which is the acute angle between the dip direction of the overall slope face along the cross-section line and the general trend of the strike-slip fault that bounds the tectonic slice. The main scope of measuring the α-angle of the tectonic slices is to evaluate its contribution to slope instabilities kinematically and to gain insight into the causes of continuous slope failures. As the α-angle decreases, the resistive effect of surface friction between the tectonic slice and the rock mass decreases, and hence, the potential instability within the tectonic slice increases. In other words, if α-angle is equal to zero, it means that the general trend of the strike-slip fault is parallel to the dip direction line of the overall slope which is the kinematically most favorable condition for sliding of the rock mass. The illustration of how to measure the α-angle is given in Figure 3.

The highest α-angles among all the tectonic slices were obtained for shear zones 1, 2, and 8 (25° $\ge$ α) which stayed stable during the operational life of the open pit (

Table 1. The α-angles for the shear zones and characteristics of shear zones.

Slice No	α-Angle (°)	Other Characteristics
Shear zone-1	25	It is in the northernmost part of the open pit.
Shear zone-2	25	These slices are narrower than the others.
Shear zone-3	10
Shear zone-4	12	The strike-slip faults form convergent lateral release surfaces for a possible mass movement.
Shear zone-5	West of the slice: 14	The western part of this zone is a potential landslide area. Hence, α: 14° was used.
	East of the slice: 23
Shear zone-6	North of the slice: 21	A large-scale landslide occurred along the strike-slip fault limiting this zone from the south on 19.12.2018. Hence, α: 2° was used for correlations with the factor of safety values.
	South of the slice: 2
Shear zone-7	North of the slice: 20	The west of this slice coincides with the western slopes of the pit that are beyond the aim of this study. Hence, α: 20° was used.
	South of the slice: 2
Shear zone-8	25	It is in the southernmost part of the open pit.

). Thus, these shear zones have the lowest potential for instability. On the other hand, the lowest α-angle was measured for the southern part of shear zone-6 in which the landslide occurred on 19 December 2018. However, the α-angle in the northern part of the same shear zone along with the strike-slip fault bounding the zone from the north is 21°. Therefore, the landslide occurred along the southern boundary of the tectonic slice, not the northern one. This suggests that the strike-slip fault limiting the zone from the south is dominating the slope instabilities. Similarly, a NE-SW (N25E/70NW) trending normal fault is dividing shear zone-5 into two pieces (Figure 1). Accordingly, the α-angle is 14° on the western piece and 23° on the eastern piece (Table 1). α: 14° means that the angle between the dip direction of benches and the strike of the strike-slip fault that is restricting the zone is 76°. Consequently, an almost rectangular shape block near the pit base is formed. Tension cracks observed in the field confirmed the potential of movement on the benches of this section (Figure 1). The higher α-angle (23°) in the eastern piece loses its significance because any slope movement at the toe of this slice will trigger the movement on the upper slope due to the undercutting effect. Accordingly, it can be stated that shear zones 5 and 6 are relatively more susceptible to sliding than the other zones. It should also be noted that shear zone-5 has a potential for slope instability not only for future pit design but also for the current conditions.

Even though shear zones 3 and 4 have relatively low α-angles (10°, 12°), they do not have critical stability conditions as shear zone 5. Because the strike-slip faults that are restricting these zones are intersecting at the base of the open pit (Figure 3). In other words, the lateral release surfaces are forming convergent boundary which blocks the sliding mass. This geometry was developed owing to the curved strike of the strike-slip fault which is restricting zone 4 from the south. As a result, the zones become narrower in the west. The potential of slope stability is therefore becoming less prominent as the rock mass within these zones is enclosed from three sides by faults and no free surface is present for a slide due to wedging effect (Figure 3).

4. Geometrical Feature of the Strike-Slip Faults: Undulation Angle

The second measured geometric feature is the undulation angle of the surface trace of the mapped strike-slip fault because it will likely contribute to the overall slope instability in the long term. Besides, they control the surface friction of the landslide mass with the release surfaces. The strikes of the foliation planes and the overall slope are nearly parallel with each other, and the foliation planes do not daylight on the overall slope. In this case, geometric features of the release surfaces of tectonic slices gain importance. If the release surfaces are undulated, the slide cannot easily occur along these surfaces even if the cohesion value of the release surfaces is almost zero. Hence, the sliding potential is promoted on the surfaces having low undulation angles which accordingly have weak attachment. In addition, as the curvature (undulation) of discontinuity increases, the friction also increases and the contact between the rock and discontinuity plane (surface trace) lengthens. In this study, the undulation angles were measured on the surface traces of nine strike-slip faults that were drawn on the 1/1000 map as illustrated in Figure 4.

In Figure 4, the III–III′ cross-section line which is constructed vertically to the strike of the overall slope was illustrated. Accordingly, the α-angle of the shear zone-5 was determined as 14°. On the other hand, the undulation angle of the surface trace of the strike-slip fault that borders the shear zone-5 from the south was determined as 5.5°.

The undulation angles of the strike-slip faults are changing between 4.1° and 8.25°. The highest undulation angle was measured on the strike-slip fault bounding the shear zone-6, while the lowest undulation angle was measured on the strike-slip fault bounding shear zone-3 (

Table 2. The undulation angles (i) of the shear zones and related cross-sections.

Shear Zone	Cross- Section	i (°)	ϕeff (°)	α Angle (°)
1,2	VII–VII′	4.50	26.50	25
3	I–I′	4.10	26.10	10
4	II–II′	6.00	28.00	12
5	III–III′	5.50	27.50	14
6	IV–IV′	8.25	30.25	2
7	V–V′	5.71	27.71	20
8	VI–VI′	5.00	27.00	25

).

5. Shear Strength Parameters of the Shear Surfaces in Orthogneisses

A series of shear tests were performed on the HQ (63.5 mm) core specimens having a natural crack approximately perpendicular to the foliation planes that are extracted from the drill hole by following the ASTM D5607-16 [37] standard. Plaster was used to encapsulate the upper and lower half of the material in the portable shear box device. Five normal stresses (50, 70, 100, 200, 350 kPa) were manually exerted at an approximate strain rate of 1 mm/min and kept constant during shearing. These surfaces are natural surfaces that are highly weathered and covered by a thin clay band. As the surface is pre-sheared, it is slickensided and smooth and hence no dilation has occurred. Consequently, shear stress-normal stress relations were derived (Figure 5). The internal friction angle (ϕ) and cohesion value (c) of these surfaces were obtained as 22° and 20 kPa, respectively.

The ϕ value obtained from the laboratory tests (ϕ_lab) was used to calculate the effective internal friction angle (ϕ_eff) as given below;

ϕ_{eff =} ϕ_lab + i

(1)

where i is the undulation angle that is used in the equation for long-term stability evaluation. ϕ_eff value for each main shear zone is presented in Table 2. The ϕ_eff values are changing between 26.1° and 30.25° for the seven geotechnical cross-sections. These values were then implemented into the numerical models while the shear strength parameters of foliation planes were estimated by the back analysis of the old landslide.

6. Rock Mass Characterization of the Orthogneisses

Two different discontinuity types have been distinguished in the orthogneiss rock mass, namely; foliation surfaces and tectonic joints. Joint sets 1 and 2 (J₁ and J₂) are “wide-spaced”, the joint set-4 (J₄) is “very wide-spaced” and the joint set-3 (J₃) is “moderately wide-spaced” according to the classification method suggested by ISRM [38]. Joint set-4 is developed nearly parallel to the normal fault-I. In addition, tension cracks are formed along J₄ in the field (Figure 6). Foliations are nearly planar or slightly undulated, moderately weathered, having an aperture of 1–5 mm with soft fillings in the sections close to the slope crest. They are also moderately wide-spaced discontinuities with high persistence.

The disturbed and highly weathered upper layers of orthogneiss that overlay the moderately weathered (MW) orthogneiss unit and the ore material were modeled in numerical simulations. The rock mass properties of MW orthogneiss and ore material were retrieved from Kadakci [31] and Karagöz and Koca [39] while the properties of the upper disturbed material in which the landslide occurred were estimated by the back-analysis of the old landslide. In this context, ranges of Hoek-Brown intact rock parameter (m_i), and GSI values were considered. Afterward, the equivalent Mohr-Coulomb parameters (c_m, ϕ_m) and modulus of deformation of the rock mass (E_m) were calculated according to the equations proposed by Hoek [40] (Figure 7). In Figure 7, the relationships between the selected (GSI, m_i) and the calculated (c_m, ϕ_m, and E_m) parameters were determined.

7. Numerical Analyses

Stress-deformation analysis or namely numerical analysis assumes deformation, stresses, and pore pressures in conjunction with the factor of safety. In the study area, multiple sets of discontinuities and close spacing of discontinuities enable the rock mass to show an equivalent near-isotropic behavior. However, failure may start from a weak discontinuity plane and then pass through a very weak rock mass [41]. Because the shear strength of rock material and the discontinuity plane may both break down if the low-strength rock material itself is accompanying the poor rock mass because of intense jointing. Numerical studies to reveal the effect of discontinuity properties such as orientation, spacing, and persistence in rock medium have been carried out considering the Mohr-Coulomb and Barton-Bandis failure criteria by several researchers using 2D finite element method (FEM) [42,43,44,45]. Hence, in this study, a combined continuum-interface method was used. This method relies on the explicit simulation of interfaces between joint walls through the joint elements introduced into the FE model running with the shear strength reduction (SSR) technique [46]. These interfaces are acting like detachment planes with tensile failure during the numerical calculation process. The shear strength of these joint elements can be readily represented by the Mohr-Coulomb criterion by using the c, ϕ, and tensile strength (σ_t). The strength reduction technique which enables calculating the strength reduction factor (SRF) which is eventually the same as the factor of safety value was used. The seven geotechnical cross-sections that are perpendicular to the overall slope and passing through each shear zone have been constructed (Figure 1). Even though there is no exact standard in the literature for choosing a domain size, the cross-sections were extended behind the pit area (>50 m) and the whole model depth was designed to be approximately 1.5 times the overall slope height. The mesh configuration was set to six noded, graded triangular elements. In order to verify the mesh size, several trial and error cycles for various mesh sizes were made until the results became insensitive to the number of elements used for the analysis [47,48,49,50]. The contiguity and quality of the mesh elements were checked based on the aspect ratio and orthogonality of the elements. The number of mesh elements are varying between 10,865 and 15,990, the lowest for VII-VII cross-section and the highest for VI-VI cross-section. Zero bad elements were detected in the models. The simulated distance between consecutive nodes which affects the precision of the calculations is 0.4 m near the ground surface and increases towards deeper elevations and reaches approximately 4.2 m. The accuracy of the models has been verified by comparing them with the in-situ displacement measurements.

Although it is impossible to model every discontinuity, engineering judgment must be used to choose which structures are the most crucial for the study area [4]. Accordingly, the joint system for the overall slopes in the NE part of the open pit was modeled by using F₁ and J₁ sets. Possible slides in the NE of the open pit will be controlled by these discontinuities. Discontinuities of J₁ longitudinally intersect the overall slope in the northeastern part of the open pit, and they form the rear release surfaces of the sliding planes. On the other hand, possible slides in the SE of the open pit will be controlled by the compression joints (J₄) and the F₁. Discontinuities of the J₄ intersect the overall slopes in the SE of the open pit and they form the tension cracks on the slope crest. Joint networks for J₁, J_4, and F₁ planes and random joints for normal faults that obey the Mohr-Coulomb failure criterion were constructed in the numerical models. The orientation, spacing, persistence, and joint end properties of the joint sets were defined. The moderately weathered orthogneiss rock mass and highly weathered, disturbed orthogneiss were modeled by the Generalized Hoek-Brown criterion [50] and Mohr-Coulomb Criterion, respectively.

7.1. Back Analyses of the Old Landslide

Shear zone-6 is the most unfavorable zone for slope stability due to the geometric conditions of the structural constraints, low rock mass quality as well as being charged by surface and groundwater, and its intrinsic permeability. Slope material of 53 m height slid towards the base of the open pit (elevations between 468 and 415 m) and the overall slope angle was 28° when the landslide occurred in 2018. The failure involved 700 × 10³ m³ rock material. The displacements progressed during a heavy rainfall period that caused high water pressures. An overburden removal at the toe was operated and the displacements were subsequently increased to 10 cm/day [34].

Pre-failed slope (IV–IV′ cross-section) was constructed using the topography measured before the landslide. It was then subjected to back analysis to determine the shear strength of highly weathered, disturbed gneiss unit and discontinuities precisely and to estimate discontinuity configuration at the time of failure (Figure 8). The back analyses method involves a probabilistic approach to the parameters that are required to be obtained [51,52,53,54]. They are based on simply finding out the combination of material properties satisfying the limit equilibrium, theoretically a factor of safety (SRF) of 1.0. However, if the structural model itself is also aimed to be determined, the back analyses need to be performed manually for trial sets of models until reaching the known conditions of the landslide. Hence, in this study, the trial numerical models were iterated until achieving the best harmonizing results with the in-situ (real) displacement values (

Table 3. The back-calculated geo-mechanical inputs.

Properties	Moderately Weathered Gneiss	Ore (Albite)	Properties	Highly Weathered, Disturbed Gneiss
	Generalized Hoek-Brown Criterion			Equivalent Mohr-Coulomb Criterion
γ (kN/m3)	26	27	γ (kN/m3)	23
σci (MPa)	27.34	75	σt (MPa)	0
mb	0.974	1.941	c (MPa)	0.5
s	0.0003	0.0015	ϕ (°)	23
a	0.508	0.504	Em (MPa)	945
Em (MPa)	2621.66	7506.61	υ	0.3
υ	0.25	0.25	-	-

σ_ci: UCS value of intact rock, σ_t: tensile strength, γ: unit weight, m_b, s, a: Hoek-Brown’s constants, c: cohesion, Φ: internal friction angle, E_m: elastic modulus of rock mass, υ: Poisson’s ratio.

and

Table 4. The back-calculated discontinuity properties.

Properties	Joint		Foliation
	J1	J4	F1
Dip angle (°)	82	80	34
Spacing (m)	10	10	5
Length (m)	50	50	10
Persistence	0.85	0.85	0.8
Joint ends	Open at boundary contact
c (MPa)	0.03		0.03
ϕb (°)	Changing between 26.1° and 30.25° depending on the shear zone		30°
kn (GPa/m)	100		100
ks (GPa/m)	10		10

k_n: normal stiffness, k_s: shear stiffness.

). For this reason, the procedure is dependent on expertise, and it is not a time-effective process. Nevertheless, it is the only way to estimate the slope conditions at the time of failure. To simulate the prolonged rainfalls as a triggering mechanism, the water table level (WTL) was kept at 470 m which is the highest level observed in the pit in winter. In addition, pore pressures were automatically calculated along the discontinuities through an effective stress approach in the numerical calculations.

The optimum model that satisfies the limit equilibrium has given a combined failure mechanism similar to the in-situ observations (Figure 9). It was determined that a two-block failure surface is initiated from the tension cracks developed from J₄ at the elevations of approximately 470 m and then breaks through the very weak rock material. Afterward, the rupture surface follows the foliation planes, which is probably the contact between the very weak gneiss and the underlying relatively stronger, moderately weathered gneiss. And finally, bulges at the elevation of 418 m. However, it should also be noted that this type of failure shows some rotational character due to the failure of weak rock layers at the point where the shear stress is the greatest along the failure surface (almost near the toe). The initiation point of the failure behind the fault plane suggested the idea that the creeks provide water flow through the joint planes which generate weak structures to cause rupture from the rock mass (Figure 9). The maximum total displacement of 6.04 cm represents the beginning of the progressive stage of mass movement at which the SRF is 1.0. At that point, the displacements have been measured as 10 cm by the data recorded by the robotic total station [34]. This shows that the 2D finite element model is in good harmony with the real case. Afterward, the slope movements increased immediately to 160 cm according to field measurements which represents the progressive phase of the slope movements.

7.2. Present Slope Stability Investigations

After the rock mass failure in 2018, the limit equilibrium conditions were achieved for IV–IV′ cross-section, and the overall slope stayed stable with a slope angle of 25°. In 2021, the slope movement again reactivated in a very low degree of displacement. The movement was measured along the inclinometer (INC-7) installed in December 2019. The sequential inclinometer measurements have been carried out four times until January 2021. The normal fault-I that is inclined at 75° toward the open pit (NW) was intersected at 19 m depth from the ground surface along the inclinometer casing. At this point, a horizontal displacement of approximately 4 mm was detected (Figure 10).

Afterward, seven geotechnical cross-sections were investigated by FE calculations considering static and pseudo-static (dynamic) conditions. It is known that WTL changes seasonally and shows fluctuations. It reaches 470 m in the winter and spring seasons when the surface flow is high. In summer, the WTL drops to 455 m. The flow direction is towards NE for slices 1, 2, 3, and 4. On the other hand, the direction of water flow is towards the base of the open pit in slices 5, 6, 7, and 8 due to the hydraulic gradient. Pore water pressure was not added to the discontinuity planes for the cross-sections passing through slices 1, 2, 3, and 4. However, it was taken into account in the rest of the cross-sections. The seismic effect on the slope stability was modeled as importing a pseudo-static load termed the seismic coefficient. The seismic coefficient was assumed to be 0.08 g. The current topography was generated from the image acquired from an unmanned aerial vehicle (UAV). The geometry and the boundary conditions are shown in Figure 11, whereas the mechanical properties and the numerical parameters adopted are listed in Table 3 and Table 4.

The SRF values derived from the FE calculations for all cross-sections are given in

Table 5. SRF values of the overall slopes along the cross-sections under static and pseudo-static conditions (with water table level at 470 m, Present).

PRESENT
Cross-Section	Shear-Zone	I (°)	ϕeff (°)	Slope Angle (°)	Hslope (m)	SRF (Static)	SRF (Pseudo-static)
I–I′	3	4.10	26.10	31	150	1.21	0.89
II–II′	4	6.00	28.00	31	144	1.12	0.85
III–III′	5	5.50	27.50	30	132	1.09	0.98
IV–IV′	6	8.25	30.25	25	123	1.62	1.27
V–V′	7	5.71	27.71	21	104	1.72	1.26
VI–VI′	8	5.00	27.00	22	105	1.70	1.29
VII–VII′	1–2	4.50	26.50	31	97	1.65	1.28

. It can be referred from Table 5 and the field conditions that the open pit slopes remain stable at present without an earthquake disturbance. The lowest SRF values for static conditions were obtained for II–II′ and III–III′ cross-sections (SRF ≤ 1.15). On the other hand, the SRF values of I–I′, II–II′, and III–III′ cross-sections are obtained as lower than 1.0 (SRF ≤ 1) under pseudo-static conditions. Because the overall slope angles and slope heights are relatively higher than the other cross-sections.

The location and depth of the inclinometer casing (INC-7) are given on V–V′ cross-section (Figure 12). It is seen that the horizontal displacement (4.7 mm) detected at 19 m depth in the FE model is almost the same as that derived from the INC-7 measurement (4 mm). This finding verified the numerical solutions.

7.3. Future Slope Stability Investigations

Since further ore production is being planned, a 25 m deepening from the current base and 25 m push-back operations have been modeled (Figure 13). The seven cross-sections have been examined under static and pseudo-static (dynamic) conditions. Different hydrogeological scenarios with varying static WTLs (470, 460, 455 m) and pore water pressures acting on joints have been employed The SRF values derived from these analyses are given in

Table 6. SRF values of the overall slopes along the cross-sections under the static and pseudo-static conditions (with water table level at 470 m, future).

FUTURE
Cross-Section	Shear-Zone	i (°)	ϕeff (°)	Slope Angle (°)	Hslope (m)	SRF (Static)	SRF (Pseudo-static)	Bench Geometry
								W (m)	H (m)	α (°)
I–I′	3	4.10	26.10	26	155	1.19	0.90	5	10	33
II–II′	4	6.00	28.00	26	182	1.28	0.95	5	10	33
III–III′	5	5.50	27.50	25	175	1.25	0.89	5	10	36
IV–IV′	6	8.25	30.25	26	169	1.40	1.05	5	10	36
V–V′	7	5.71	27.71	25	131	1.32	0.94	5	10	37
VI–VI′	8	5.00	27.00	26.5	126	1.28	0.97	5	10	38
VII–VII′	1–2	4.50	26.50	30	109	1.50	1.27	5	10	37

,

Table 7. SRF values of the overall slopes along the cross-sections under the static and pseudo-static conditions (with WTL lowering to 460 m, future).

FUTURE
Cross-Section	Shear-Zone	i (°)	ϕeff (°)	Slope Angle (°)	Hslope (m)	SRF (Static)	SRF (Pseudo-static)	Bench Geometry
								W (m)	H (m)	α (°)
I–I′	3	4.10	26.10	26	155	1.25	0.96	5	10	33
II–II′	4	6.00	28.00	26	182	1.38	1.04	5	10	33
III–III′	5	5.50	27.50	25	175	1.35	0.95	5	10	36
IV–IV′	6	8.25	30.25	26	169	1.47	1.11	5	10	36
V–V′	7	5.71	27.71	25	131	1.49	1.05	5	10	37
VI–VI′	8	5.00	27.00	26.5	126	1.32	1.01	5	10	38
VII–VII′	1–2	4.50	26.50	30	109	1.63	1.29	5	10	37

and

Table 8. SRF values of the overall slopes along the cross-sections under the static and pseudo-static conditions (with WTL lowering to 455 m, future).

FUTURE
Cross-Section	Shear-Zone	i (°)	ϕeff (°)	Slope Angle (°)	Hslope (m)	SRF (Static)	SRF (Pseudo-static)	Bench Geometry
								W (m)	H (m)	α (°)
I–I′	3	4.10	26.10	26	155	1.35	0.98	5	10	33
II–II′	4	6.00	28.00	26	182	1.40	1.06	5	10	33
III–III′	5	5.50	27.50	25	175	1.37	1.02	5	10	36
IV–IV′	6	8.25	30.25	26	169	1.51	1.17	5	10	36
V–V′	7	5.71	27.71	25	131	1.56	1.11	5	10	37
VI–VI′	8	5.00	27.00	26.5	126	1.38	1.06	5	10	38
VII–VII′	1–2	4.50	26.50	30	109	1.67	1.35	5	10	37

.

Among the all cross-sections, III–III′ cross-section stays within shear zone-5 and intersects the potential landslide area close to the current base of the open pit (432–371 m). The tension cracks mapped on a 1/1000 scale within the shear zone-5 are observed as starting from the 7th bench (430 m) from the current base of the mine. As the current base is 25 deepened for further ore production, the failure may be expected to start from the 9th bench (430 m). However, Figure 14b indicates failure surface can be retrogressive due to the undercutting effect and may start from 470 m in the α:23 area.

The decrease in WTL level from 470 m to 455 m along the cross-sections caused a decrease in pore pressures. Therefore, the decrease in WTL led to an average increase in SRF values of 11.09% and 11.45% under static and dynamic conditions, respectively (Table 6 and Table 7). The greatest increments in SRF value both for static and dynamic conditions were obtained for V–V′ cross-section (18.18%, 18.09%), while the lowest increments were for VII–VII′ (11.33%, 6.29%) cross-section. On the other hand, the highest SRF differences between dynamic and static conditions as WTL is lowered were obtained for I–I′, III–III′, VII–VII′ cross-sections. It means that these cross-sections tend to be more sensitive to dynamic conditions while being water saturated. Besides, the lowering of WTL from 460 m to 455 m led to an increase in SRF values an average of 3.62% and 4.58% under static and dynamic conditions, respectively. The increase in SRF values under the dynamic conditions due to WTL lowering of 5 m (460 m to 455 m) was determined between 1.9% and 7.4%. These findings implied that the lowering of WTL from 460 m to 455 m has a more pronounced effect under dynamic conditions.

8. The Relationships between the Geometrical Features, WTL, and SRF Values

Joint water pressure is important in terms of slope stability, especially for deep rock slopes with a critical balance. As in this work, the change in WTL between 460 m and 455 m significantly affected the stabilities of cross-sections with the same geometric and seismic conditions. The slope facets are unevenly oriented due to the different shear zones present in the open pit. The concave and convex slope geometries as well as convergent and divergent shear zones are observed. Nian et al. [20] determined that the concave-shaped slope face yields higher SRF than the convex one. This finding relies on the different strain distributions. On the concave slope facet, the shear strain mainly distributes along the limbs while concentrating in the center on the convex slope facet [20]. Hence, the stress accumulation within the pit principally concentrates in the shear zones 4, 5, and 6. On the other hand, the concave slope facet accelerates water flow inward which means water accumulates in the center of the slope [55]. In addition, the normal faults dipping towards the pit with different trends also provide surface water flow into the pit. That is why several cross-sections investigated in this study are highly affected by WTL lowering due to having a concave slope facet.

Except for the slices including VII–VII′ cross-section, all the other cross-sections yielded critical balance (SRF~1.0) in terms of the overall slope stability under dynamic conditions. Nevertheless, if WTL is kept at an elevation of 455 m, the proposed slope angles can provide stable conditions (Table 8). Again, except for VII–VII′ cross-section, the critical conditions are achieved for the overall slope angles between 25° and 26.5° for all other geotechnical cross-sections. The proposed overall slope angles are a bit far from the economic limit for mining operations. Accordingly, the overburden stripping ratio will significantly decrease which directly increases the ore-production cost. It was also determined by the detailed fracture mapping that the surface water recharges the groundwater in rainy seasons for the shear zones 5, 6, and 7. Besides, the water recharging due to intense drainage patterns is more pronounced in the shear zone-6. The embedded stream bed in the south of this slice that is aligned parallel to the strike-slip fault carries all stream water and works as a drainage path for discharging into the open pit. Consequently, the landslide occurred on 19 December 2018 along this embedded stream bed and disrupted the mining operations.

As it was previously mentioned, WTL lowering affected the SRF values principally under dynamic conditions. Therefore, the decrease in WTL from 470 m to 455 m was examined graphically with the H_slope versus SRF curve (Figure 15). Ranges of variation in SRF values for I–I′ (150 m), II–II′ (182 m), III–III′ (175 m), IV–IV′ (169 m), V–V′ (131 m), VI–VI′ (126 m), and VII–VII′ (109 m) were obtained as 0.08, 0.11, 0.13, 0.12, 0.17, 0.09, and 0.08, respectively. It was understood that the cross-sections having concave slope facets are more sensitive to changes in WTL than the other sections. Similar results can be declared for static conditions.

One of the geometrical features proposed in this study that affects the slope stability condition in the field is the α-angle. α-angles were measured on the map with a 1/1000 scale because it is not possible to accurately measure these angles in the field. α-angle can be served as an effective parameter for assessing the sliding potential of tectonic slices kinematically. It is observed in the open pit that as α-angle increases, the ability to resist sliding also increases. This argument was verified by evaluating the relationship between the α-angle and SRF values with varying WTLs (Figure 16). The yellow-ellipse area on the graph given in Figure 16 represents the IV–IV′ cross-section on which the landslide occurred previously. It can be referred from Figure 16 that as the α-angle increases, the SRF values tend to increase. Even though IV–IV′ cross-section is in the α:2° area, after the landslide that moved the 700 × 10³ m³ material downward, the overall slope angle decreased to 25° from 28°, and this cross-section remains stable at 25° in the present condition. Therefore, SRF values are obtained higher than expected (Figure 16). This graphical finding also supports that this cross-section has reached the limit-equilibrium condition. The SRF values of the other cross-sections showed a good correlation with α-angles. A strong relationship is observed both for static and dynamic conditions.

Another geometrical feature proposed in this study that affects the slope stability condition is the undulation angle (i). Therefore, the effect of undulation angle on the slip potential of the tectonic slices was investigated by implementing these values into FE calculations with the ϕ_eff = ϕ_lab + i equation. i-angle is the undulation angle measured along with the surface trace of the strike-slip fault. The changes in SRF = a × e^{−b (H)} curves were examined for the cases if the ϕ_eff value is increased by 1°, 2°, and 2.5°, (Figure 17). These investigations were performed separately for the WTLs of 470 m and 455 m. It was determined that as the undulation angles of the strike-slip faults increase, SRF values of the overall slopes also increase. Again, the SRF values of the IV–IV′ cross-section are discordant with the general trend (Figure 17). That is why the correlation coefficients referred to fair relationships rather than strong relationships for future scenarios. If the slope heights and slope angles of each cross-section are kept constant and the ϕ_eff value of the strike-slip faults is increased gradually by 2.5°, the changes in the SRF values for the water table levels of 470 m and 455 m were determined to be similar (0.05–0.11 and 0.04–0.11, respectively). However, Figure 17 shows that the range of these changes decreases as the slope height decreases. This finding suggests that the positive effect of ϕ_eff value on the slope stability is decreasing as the slope height decreases. Nevertheless, these increments in SRF values are still crucial for the overall slopes that have critical stability.

The lowest undulation angle was measured for the shear zone-3. Accordingly, the SRF values obtained for the I–I′ cross-section under static and dynamic conditions for the WTL of 455 m were found to be less than the other cross-sections. On the other hand, the SRF values obtained for the II–II′, III–III′, VI–VI′ that is within the shear zones 4, 5, and 8 are approximately 1.0 under dynamic conditions. Their α-angles are between 2° and 25°, while i-angles are between 5° and 6°. It is likely therefore that not only the α-angles and i-angles are controlling the stability but also the slope height, slope angle and pore pressures govern the slope instabilities.

To sum up, the variations in hydraulic gradient, intense drainage pattern due to shearing, and concave slope facet control the influence of water on slope stability. Besides, the undulation angle of dominant discontinuity planes (strike-slip faults in this study) directly affects the shear strength and hence the slope stability while the α-angles are contributing block formation. These blocks can easily detach from the mass in the case of lower α-angles. Overall slope angle and slope height also influence the amount of shear accumulation at the toe within the time. Thereafter, the shear accumulation at the toe triggers the slope instabilities.

9. Conclusions

In the study area, mass displacements with different scales have been continuing due to varying environmental conditions during ore production since 2010. More recently, a large-scale landslide occurred on the southeastern slopes of the open pit on 19 December 2018. Even though remedial measures have been operated, the in-situ monitored displacements began to increase in 2021. In addition, the extraction of ore from deeper levels of the pit from the current base is planned which will most probably cause unstable conditions. For these reasons, it is required to determine the root causes of slope instabilities.

Back analysis of the pre-failed slope provided a realistic engineering geological model for the open pit mine. It was determined that the old landslide is a two-block, combined failure of the heavily jointed rock mass releasing from the dominant discontinuity planes. It can be argued that the stability of the slopes is affected to a certain extent (up to 18.2%) by the undulation angles (i) of the lateral release surfaces and accordingly the ϕ_eff value. This is the first time that a pit-scale undulation angle was used to investigate the stability of open pit slopes. This finding is very significant for deep rock slopes in a critical balance. Therefore, undulation angles of the surface traces of the pit-scale (≥270 m) structural features affecting the slope stability should be measured on the map to achieve a realistic approach. Because the boundaries of these tectonic slices are forming lateral release surfaces which is one of the main causative factors of landslides. It was also determined that the effect of joint water pressure on the stability analyses is crucial as well as the increase in ϕ_eff. On the other hand, the α-angle proposed in this study is also strongly related to the stability of the slope. This argument was verified by the correlations of α-angle with SRF values. The feasibility and effectiveness of the proposed parameters were validated by in-situ measured displacements and by the interactions with other environmental conditions such as WTL. Taken together, these findings suggest that a combined effect of several parameters controls the stability of open pit slopes. Therefore, it is recommended to rely on the geometrical features of the shear zones in open pit slope design in tectonically disturbed, highly fractured zones for sustainable mining. The parameters proposed in this study can also be extended to areas of similar geological conditions. As the number of trials by 3D numerical modeling increases, it can provide more definitive evidence for the contribution of the undulation angle of pit-scale lateral release surfaces to the slope stability evaluation.