







































Georgian Geographical Journal 

 

Assessment of the stability of some 

landslide bodies in Georgia using a 

numerical-analytical approach 
Levan Japaridze1* , Nikoloz Chikhradze1,2 , 

Tamaz Akhvlediani1 , Tamaz Gobedzhishvili1 , 
1 LEPL Grigol Tsulukidze mining institute, Tbilisi, Georgia 
2 

Georgian Technical University, Tbilisi, Georgia
 

* Corresponding author: levanjaparidze@yahoo.com 

 

 

 

 

Introduction 

The history of calculating the stability of landslide-prone slopes dates back 250 years. To date, more 

than 5000 articles or monographs have been published. Most of them are simplified methods for 

calculating the stability of landslide bodies of assequental and consequental (fig.1, a) type, i.e. without 

stratification and sliding on stratification planes. In their calculation schemes, the landslide body is 

divided into vertical blocks, the stability coefficients of the individual block and then of the entire 

landslide body are determined by solving systems of equations with many unknowns under simplified 

assumptions. The weaknesses of this type of methods were noted by the authors themselves (Janbu, 

Bishop and others). Their critical analysis is also given, for example, in (Stability Modelling with 

SLOPE/W 2012). These shortcomings remain in a number of recently created computer programs (for 

example, "Rocscience" Slide 2, Slide 3, Slide 6 and others), which are specially designed for landslides 

and do not consider the stress-deformed state analysis of the landslide body. 

In contrast to asequential and consequental, the first papers on the stability of block-flexure (Fig.1, b) 

landslides have appeared in recent years, and the first steps are still being taken to develop their stability 

calculation methods. 

 
Figure 1. Schematic images of a) consequental and b) block-flexure landslides 

Estimating the stability of block-flexure landslide bodies is much more difficult and calculation 

methods are less compared to consequential ones. As early as 1976, R. Goodman and J. Bray (Goodman 

& Bray1976) wrote that "our knowledge of this type of phenomena is still in its infancy". They compiled 

possible calculation schemes and identified the forms of collapse: flexural toppling, block toppling, 

block-flexural toppling (Goodman & Bray 1976), whose corresponding calculation methods are still in 

the study and development stage (Amini et. al., 2008). 

Georgian Geographical Journal, 2024, 4(2) 14-21 

© The Author(s) 2024 

 
This article is an open access article distributed under 

the terms and conditions of the Creative Commons 

Attribution (CC BY) license (https:// 

creativecommons.org/licenses/by/ 4.0/). 

DOI: 

https://journals.4science.ge/index.php/GGJ 

Abstract 

Numerical-analytical methods of stability of consequential, block-flexure and 

creep-plastic landslide-prone slopes are being developed at the G.Tsulukidze 

Mining Institute. These methods have been used: to determine the degree of 

stability of the landslide-prone blocks remaining after the 2015 landslides in 

Tskneti-Samadlo, to calculate the stability of the landslide bodies of M. 

Machavariani and Sheshelidze Streets in Tbilisi, and then - to develop their 

strengthening projects together with the specialists of the German-Bavarian 

Bureau of Engineering Geology and the construction company CRP. 

Keywords: numerical-analytical method; stability of consequential block-

flexure and creep-plastic slopes. 

Citation: Japaridze, L.; Chikhradze, N.; 

Akhvlediani, T.; Gobedzhishvili, T. 
Assessment of the stability of some 
landslide bodies in Georgia using a 

numerical-analytical approach. 
Georgian Geographical Journal 2024, 

4(2). 14-21 

https://doi.org/10.52340/ggj.2024.04.02.02 



Japaridze et al. 2024 4(2) 

15 
 

In the vast majority of the methods developed for the study of landslide events, the stability of the 

slope is evaluated under the conditions of limit equilibrium, Coulomb-Moir and similar strength 

theories, and the influence of the time factor can only be reflected indirectly, according to the studies 

of recent years, through the limit values of the rock adhesion and the internal friction coefficient, i.e. 

the strength indicators, at the expense of degradation. This was called the "Shear Strength Reduction 

Factor" and its use in conjunction with the "Finite Element Technique method new era in the theory of 

landslide stability calculation (Rocscience, 2004). 

Consideration of the "shear strength reduction factor" is very important when calculating the stability 

of block-type slide bodies that maintain a completely stable state for a long time and suddenly lose their 

stable equilibrium as a result of extreme meteorological challenge or other effects. Until then, long-

term monitoring of such landslide-prone slopes may not show any significant rock movement. 

In contrast, creep plastic deformation is more rheological in nature. Depending on the intensity of the 

active forces acting in a particular scheme and the mechanical parameters of the massif, creep 

deformation can occur very rapidly or very slowly, in periods ranging from a few hours to a hundred 

years. It is not only the strength (stress, angle of internal friction) but also the deformation characteristics 

(modulus of elasticity, creep parameters, etc.) of the rock that should be most characteristic. These can 

be determined by laboratory testing of samples of geological material and/or by processing images 

obtained from preliminary monitoring. "Natural analysis method (back analysis). The means for its full 

use will be created after the creation and improvement of the apparatus for the direct analysis of self-

creeping plastic landslides. In all these directions, the development of methods for calculating the 

stability of landslide bodies and anchored constructions is a scientific and technical problem of 

international level. 

Methods and Materials 

The numerical-analytical version of the method for assessing and calculating the stability of landslide 

bodies, both unstratified and consequential, was developed at the G. Tsulukidze Mining Institute. This 

was done by combining a stress analysis apparatus and computer programs, in particular Phase2, 

developed by Rocscience not for landslide problems. Unlike existing methods, this version makes it 

possible to simulate design schemes of various forms of possible collapse of landslide bodies, based on 

an analysis of the stress-strain state of both an unstratified and a consequent structure; the necessary 

measures to ensure its stability under normal and extreme conditions under the influence of gravitational 

and seismic forces must be calculated. 

 
Figure 2. The topographic plan, top view and a section of remained rock block - “Southern Cliff”. 

The computerised numerical-analytical method for determining the stability of landslide-prone blocks 

was used to ensure the safety of the remaining blocks after the 2015 landslide events in Tskneti-Samadlo 

and is described in (Japaridze et. al. 2020). Situation images of one example of them are shown in Fig. 

2. Using the calculation diagram and the Rocscience Phase2 computer program, horizontal (XX), 

vertical (YY) and shear (XY) stresses caused by gravitational and seismic forces were determined. The 

corresponding polychromes of vertical and shear stresses are shown in Fig. 3. 



Japaridze et al. 2024 4(2) 

16 
 

 
Figure 3. Finite element model of the "upper cliff" of the 2015 large landslide of Tskneti-Samadlo 

The values of engineering geological elements, cohesion (c) and angle of internal friction (φ) for the 

sliding surfaces of the blocks were determined by preliminary exploration studies (Table 1). 

Table 1. The values of engineering geological elements for the sliding surfaces 

Type of discontinuity 

Sliding plane 

Filling Cohesion c’ 

[MN/m2] 

Friction φ 

[o] 

Slickenside (normal fault) Locally fault gouge 0 15-20 

Join (rough, medium transection) Locally clay/silt (0.005) 28 

Bedding plane in claystone - 0.05 22.5 

Bedding plane in strongly 

weathered claystone/clay 

- 0.01 20 

Using Rocscience Phase2, the values of horizontal (XX), vertical (YY) and tangential (XY) stresses 

at 30 points of the potential slip planes were determined and entered into Excel Table 2. At the same 

points on the slip planes with an angle of θ, the normal N and shear T stresses are calculated according 

to the following relationships: 

𝑁 = 0.5(𝑋𝑋 + 𝑌𝑌) + 0.5(𝑋𝑋 − 𝑌𝑌)𝑐𝑜𝑠2𝜃 + 𝑋𝑌𝑠𝑖𝑛2𝜃                        (1) 

 𝑇 = 0.5(𝑌𝑌 − 𝑋𝑋) 𝑠𝑖𝑛2𝜃 +𝑋𝑌 𝑐𝑜𝑠 2𝜃                        (2)  

where θ is the angle from the positive direction of the x axis to the plane on which the stresses (1), (2) 

are determined. 

By inserting them and the rock cohesion C and the coefficient of internal friction  on the sliding 

surface into the Coulomb-Maure limit equilibrium formula (3), stability is estimated at each (30) point 

and then over the entire area of the slope under normal and extreme conditions. 

[𝑇] = 𝐶 + 𝑁𝑡𝑔𝜑 ≥ 𝑇                       (3) 

The values of the holding stresses [T] are entered into the Excel 2 table and the safety factors are 

calculated at the points of the potential slip surface (Fig. 3). For the given natural conditions C = 50 

kPa, φ = 22.50, θ = 190, the total safety factor is SF = 6,07 and stability is maintained by a large margin. 

Table 2. Excel table for determining safety factors at potential sliding surface points 

№ XX kPa YY kPa XY kPa N kPa T kPa [T] kPa SF 

1 1.3 95.4 12.1 19 4.9 15 2.97 

2 0.9 92.5 10.8 17 5.0 15 2.90 

3 0.6 86.9 9.3 15 4.9 14 2.93 

4 1.0 76.9 9.6 15 4.0 14 3.56 

5 5.8 66.9 20.7 25 0.6 15 23.93 

- - - - - - - - 

26 1.3 95.4 12.1 19 4.9 15 2.97 

27 0.9 92.5 10.8 17 5.0 15 2.90 

28 0.6 86.9 9.3 15 4.9 14 2.93 

29 1.0 76.9 9.6 15 4.0 14 3.56 

30 5.8 66.9 20.7 25 0.6 15 23.93 



Japaridze et al. 2024 4(2) 

17 
 

Due to the possible reduction of the rock strength characteristics in extreme conditions, the results of 

the calculation for the reduced rock shear strength indicators obtained as a result of additional 

investigations by the specialists of the Bavarian State Office for Engineering Geology are given in Table 

3, from which it can be seen that with zero cohesion (c=0) and the coefficient of internal friction =150. 

In this case, the safety factor is less than one (SF=0.88), i.e. the slope will not be stable. Appropriate 

geotechnical measures were therefore recommended. 

Table 3. Safety factors due to a possible reduction in the strength characteristics of rocks 

Cohesion C, kPa Internal friction  φ, (0) Safety Factor SF 

50 22.5 6.07 

10 20 2.13 

0 15 0.88 

The proposed numerical-analytical approach for calculating the stability of a block-flexure landslide-

prone slope, whose finite element model is shown in Fig. 4, reflects the block-flexural toppling 

mechanism much more accurately than the current kinematic schemes of Goodman and Bray (Goodman 

& Bray 1976). Analytical and numerical methods have been developed. It differs substantially from the 

special computer program "RocTopple 1.0" recently developed by the company "Rockscience", which 

literally repeats the idealised calculation scheme of Goodman and Bray with very important 

simplifications and assumptions about the stiffness of the blocks, the equal thickness of the layers and 

their stability shapes. As an example of non-sequential landslides, we considered the landslide-prone 

areas on the slopes of the orographically left bank of the Vere River on the Chabua Amirejibi road. The 

stability coefficient of one of the slopes against gravitational and seismic loads was determined 

(SF=16.5). 

 
Figure 4. Finite element model of block-flexure slope on Chabua Amirejibi Street  

In a large number of mechanical models describing landslide phenomena, rock masses are considered 

from the point of view of solid body mechanics. The stability of slopes is assessed under the conditions 

of limit equilibrium, the Coulomb-Moire theory and similar strength theories. Research in recent years 

has shown that the influence of the time factor can only be reflected indirectly in these methods by 

reducing the limit values for rock adhesion and the coefficient of internal friction, i.e. strength 

indicators. This has been termed the "shear strength reduction factor" and its use in combination with 

the finite element method has been described as "a new era in landslide stability analysis" (Stability 

Modelling with SLOPE/W, 2012). 

Consideration of the 'sliding strength reduction factor' is very important when calculating the stability 

of block-type consequential and non-stratified landslide bodies, which maintain a stable state for a long 

time and then suddenly lose equilibrium as a result of extreme meteorological or other effects. Until 

then, long-term mini-tours of this type of landslide-prone slope may show no significant rock 

movement. 

Creep-plastic deformations occur practically everywhere in rock massifs, and depending on the 

intensity of the active shear stresses acting in a particular scheme and the mechanical parameters of the 

massif, creep deformations can develop very rapidly or extremely slowly, in periods ranging from a few 

hours to hundreds of years. Active compressive stresses are mainly those caused by the gravitational 

forces acting in the rock mass. As far as the mechanical parameters of the massif are concerned, the 

most characteristic are not only the strength (density, angle of internal friction), but also the deformation 

characteristics of the rocks: modulus of elasticity, creep parameters. These can be determined by 

laboratory analysis of samples of geological material taken from the likely landslide-prone slope and/or 

by processing the image obtained by preliminary monitoring. using the "natural-analytical", i.e. reverse 

analysis (back analyse) method (Japaridze,1987). It should be noted that the means for the full use of 

the latter will be created in the wake of the creation and improvement of the direct analysis apparatus 

for creep-plastic landslides, which is still a current problem. 

69.15

-308.01

-265.05

-235.78

-210.46

-187.58

-156.67

-156.92

-140.87

-135.04

-122.36

-109.88

-105.55

-104.09

-112.54

-97.86

-88.04

-82.25

-77.17

-70.85

-71.92

-57.64

-55.35

-46.83

-37.18

-33.17

-29.09

  0.17

User Data

Tau XY

-425.00

-400.00

-375.00

-350.00

-325.00

-300.00

-275.00

-250.00

-225.00

-200.00

-175.00

-150.00

-125.00

-100.00

 -75.00

 -50.00

 -25.00

   0.00

  25.00

  50.00

  75.00

40
30

20
10

0

0 10 20 30 40 50 60 70 80 90 100



Japaridze et al. 2024 4(2) 

18 
 

It is difficult to assess the stability of creep-plastic landslide bodies in terms of completeness, accuracy 

and practical possibility of taking into account the factors, and there are practically no calculation 

methods. Therefore, the modern stage of development of mining massif rheology can be considered 

only the beginning, when the primary stock of information, ideas, methods and results has been 

accumulated, but the main researches and practical application of their results are still ahead. 

In this case, a study on the stability assessment of creep-plastic landslide slopes was considered in 

order to develop a calculation method which, unlike the method for calculating the stability of a block-

type landslide body, should take into account the elastic-plastic and rheological deformations occurring 

continuously in the massif containing metamorphic clay rocks and the time factor. This problem is 

important and relevant both theoretically and practically in general, and also because there are 

geological formations in many places in Georgia where we have "creep-plastic" landslide-prone places. 

Local or large-scale landslides of this type are to be expected in areas with hilly terrain (Ministry of 

Environmental Protection and Agriculture of Georgia 2019) if, under the conditions of intensification 

of infrastructure construction, special attention is not paid to their prevention even at the early stage of 

design. 

At present it is possible to determine the functions representing the displacement of the points of the 

array, to assess the stress-strain state of geotechnical objects taking into account creep, if the analytical 

solution of the task within the theory of elasticity is known. Such tasks include underground structures, 

namely tunnels, which can be mathematically modelled, for example, using the methods of complex 

potentials of N. Muskhelishvili. However, it is very difficult to do the same when assessing the stability 

of slopes prone to landslides. 

The use of finite element methods in a limit equilibrium framework to analyse the stability of 

geotechnical structures is a significant step forward and has many advantages: 

− It is no longer necessary to make assumptions about internal stresses (between prisms); 

− The stability factor is determined after the stresses have been calculated, so there are no 

iterative convergence problems; 

− The problem of mobility compatibility is solved; 

− Calculated stresses are closer to reality; 

− Stress concentrations are considered indirectly in the stability analysis; 

− Dynamic stresses resulting from an earthquake can be considered directly under stability 

conditions. 

Thus, the finite element approach overcomes many of the limitations inherent in limit equilibrium 

methods. The tools needed to perform stability analysis of geotechnical structures based on stresses 

calculated by the finite element method are available today, but they also have some weaknesses, and 

unforeseen problems may arise in the future. Such problems are currently the subject of research and 

will be solved over time as the method is increasingly used in geotechnical practice (Stability Modelling 

with SLOPE/W, 2012). 

Due to the great practical importance of the task and the theoretical difficulties of a strict solution, the 

modelling of the creep-plastic landslide body can also be considered as a so-called "variable module" 

(Nicholson et. al. 1990; Technical Manual for Design and Construction of Road Tunnels, 2009) 

The results of the research were used in the creep-plastic landslide reinforcement project of Tbilisi's 

Vashlijvari district, prepared by Caucasus Road Project LLC and the Bavarian Office of Engineering 

Geology, whose plan and section C-C as of the beginning of March 2021 is shown in Fig. 5. 

 
Figure 5. Top view and plan of "Vashlijvari Landslide" with observation points 

To detect morphological changes in the landslide area, scientists from the University of Vienna carried 

out observations and satellite measurements using the Sentinel satellite's Interferometric Synthetic 



Japaridze et al. 2024 4(2) 

19 
 

Aperture Radar (InSAR). Observation points were placed on the landslide body and beyond (Fig. 5). 

The displacements of the geodetic markers in the direction of the landslide were recorded for 100 days. 

 

 

 
Figure 6. Movements of observational geodetic markers in the direction of the landslide 

The displacements caused by the gravitational forces of the observational geodetic markers shown 

in Fig. 6 make it possible to determine the creep parameters of the entire landslide body massif using 

the proposed natural-analytical method, the so-called use of time-varying modules. 

 

 
Figure 7. Polychromy of stresses XX and XY acting in the CC plane when anchors are placed at two levels 

 
Figure 8. Project of anchor reinforcement of the Vashlijvari landslide slope and its fragment 

These data, obtained by the natural analytical method carried out on the landslide body of M. 

Machavariani Street, were used for the calculation of the stability of the landslide body, for the 

preparation of the strengthening project (Fig. 7, Fig. 8) and its implementation. 

Results 

1. The original numerical-analytical version of the method for assessing and calculating the stability 

of landslide bodies, both unstratified and consequential, was developed. Unlike existing methods, this 

version makes it possible to simulate design schemes of various forms of possible collapse of landslide 

bodies, based on an analysis of the stress-strain state under the influence of gravitational and seismic 

forces; 

2. Proposed numerical-analytical approach for calculating the stability of a block-flexure landslide-

prone slope reflects the toppling mechanism much more accurately than the current kinematic schemes. 

As an example of non-sequential landslides, landslide-prone areas on the slopes of the orographically 

left bank of the Vere River were considered. 



Japaridze et al. 2024 4(2) 

20 
 

3. Morphological changes and displacements of geodetic marks in the direction of a creeping-plastic 

landslide body in the Vashlijvari region of Tbilisi were established together with scientists from the 

University of Vienna using the interferometric synthetic aperture radar (InSAR) of the Sentinel satellite. 

These data were used to calculate the stability of the landslide body using the natural-analytical method, 

to prepare a project for strengthening the slope and its implementation. 

Conclusion 

When checking the condition of stability of landslide-prone slopes and calculating their Safety Factor 

(SF), the rock mass should be characterised by its long-term shear strength and full deformation 

expected in all parts of the sliding plane. In cases where the laboratory study of the engineering-

geological material of the object as a result of practically instantaneous tests gives indicators of the peak 

strength of the rocks, the stability condition cannot be lower than the ratio of the immediate and long-

term strength; 

It is well known that the rock mass passes through different stages of deformation before collapsing 

and undergoes significant plastic deformations, the magnitude of which depends on lithology, density, 

moisture, liquefaction during seismic impact and many other indicators. All the rocks that make up the 

stratified massif and that participate in the landslide event are characterised by different deformation 

properties. Therefore, even if they undergo almost the same deformation at each moment, the stresses 

acting on them are different. In the less plastic layer, the shear deformation stresses reach the limit, 

while in the more plastic layer they are still far from the limit. The magnitude of the resistance to slope 

displacement is not equal to the sum of the resistances of the individual layers at any stage of 

deformation;  

The value of the creep-plastic landslide stability coefficient depends not only on the Coulomb-Mohr 

strength characteristics, but also on the deformation parameters of the rocks in general, including the 

creep characteristics; 

The treatment of displacements caused by the gravitational forces of observational geodetic marks by 

the proposed "natural-analytical method" compared to the cylindrical dilatometric method indicated in 

the normative documents of the Eurocodes, provides a better opportunity to determine the equivalent 

creep parameters of the entire landslide body array. 

Due to the great practical importance of the task and the theoretical difficulties of a strict solution, the 

possibility of modelling a creep-plastic landslide body should be studied by the so-called Using the idea 

of a "variable module". In this case, it is necessary to consider the proposal and development of the 

"Shear Modulus Reduction" factor similar to the already established "Shear Strength Reduction" factor; 

The value of the creep-plastic landslide stability coefficient depends not only on the Coulomb-Mohr 

strength characteristics, but also on the deformation parameters of the rocks in general, including the 

creep characteristics; 

The treatment of displacements caused by the gravitational forces of observational geodetic marks by 

the proposed "natural-analytical method", compared to the cylindrical dilatometric method indicated in 

the normative documents of the Eurocodes, provides a better opportunity to determine the equivalent 

creep parameters of the entire landslide body array. 

Due to the great practical importance of the task and the theoretical difficulties of a strict solution, the 

possibility of modelling a creep-plastic landslide body using the so-called "variable modulus" idea 

should be studied. In this case, it is necessary to consider the proposal and development of the "shear 

modulus reduction" factor, similar to the already established "shear strength reduction" factor; 

geological investigation of a specific landslide-prone site, together with other rock parameters, "natural-

analytical" and/or laboratory methods should determine rock creep parameters, or in extreme cases - 

peak and residual deformation modulus. The quantities are as required, for example, in the "Material 

Properties" table of the Phase 2.7 program. It should also be noted that the magnitude of the "time-

varying modulus" of deformation depends on the stress-strain state of a particular rock mass, i.e. on the 

overall calculation scheme of the geotechnical situation. Similarly, the creep characteristics will not be 

invariant to stress, since their curves reflect the evolution of both creep and creep-plastic deformation 

over time. 

The proposed natural-analytical method has been used: to determine the degree of stability of the 

landslide-prone blocks remaining after the 2015 landslides in Tskneti-Samadlo; St. to calculate the 

stability of the landslide bodies of M. Machavariani and Sheshelidze Streets in Tbilisi, and then - to 

develop their strengthening projects together with the specialists of the German-Bavarian Bureau of 

Engineering Geology and the construction company CRP. 



Japaridze et al. 2024 4(2) 

21 
 

Competing interests 

The authors declare that they have no competing interests. 

Authors’ contribution 

L. J. conceived of the presented idea and performed the numerical-analytical calculations, N. C. 

conceived of the presented idea, T. A. provided critical analysis, T. G. provided critical analysis and 

performed the analytic calculations. 

ORCID iD 

Levan Japaridze https://orcid.org/0000-0002-1363-0883 

Nikoloz Chikhradze https://orcid.org/0000-0001-5369-5123 

Reference 

Amini, M., Majdi, A., & Aydan, Ö. (2008). Stability analysis and the stabilisation of flexural toppling failure. 

Rock Mechanics and Rock Engineering, 42(5), 751–782. https://doi.org/10.1007/s00603-008-0020-2 

Goodman, R.E., & Bray, J.W. (1977). Toppling of Rock Slopes. 

Japaridze, L., Neumann, P., & Trapaidze, P. (2020). A Case Study of the Slope Stability after Large Landslide in 

the 2015 Flood in Tbilisi. Tbilisi. 

Japaridze, L. (1987). Analytical apparatus for determining rheological parameters by measuring the displacement 

of the tunnel wall of a circular section driven in an array of elastic-viscous plastic rocks. Bulletin of the 

National Academy of Sciences, Tbilisi. 

Ministry of Environmental Protection and Agriculture of Georgia, National Environmental Agency. Department 

of Geology. Engineering-geodynamic conditions of the territory of Tbilisi and assessment of geological 

hazards, (2019). Tbilisi, PDF. 

Nicholson, G.A. & Bieniawski, Z.T., (1990). A nonlinear deformation modulus based on rock mass classification. 

Int. J. Min. & Geological Engineering 8, 181-202. 

Rocscience, (2004). A New Era in Slope Stability Analysis: Shear Strength Reduction Finite Element Technique. 

Strength reduction. https://static.rocscience.cloud 

Stability Modelling with SLOPE/W. An Engineering Methodology, (2012). GEO-SLOPE International 

Ltd.Calgary, Alberta, http://www.geo-slope.com. 

Technical Manual for Design and Construction of Road Tunnels - Civil Elements, (2009). US Department of 

Transportation Federal Highway Administration, Publication No. FHWA-NHI-10-034. December. 

https://doi.org/10.1007/s00603-008-0020-2

