=Paper=
{{Paper
|id=Vol-2608/paper43
|storemode=property
|title=Parallel simulation of rock deformation processes
|pdfUrl=https://ceur-ws.org/Vol-2608/paper43.pdf
|volume=Vol-2608
|authors=Yevhen Tiurin,Antonina Andrieieva
|dblpUrl=https://dblp.org/rec/conf/cmis/TiurinA20
}}
==Parallel simulation of rock deformation processes==
https://ceur-ws.org/Vol-2608/paper43.pdf
Parallel simulation of rock deformation processes
Yevhen Tiurin1[0000-0001-8362-9223], Antonina Andrieieva1[0000-0002-0361-5436]
1
Donetsk National Technical University, Shybankova sq., 2, Pokrovsk, 85300, Ukraine
yevhen.tiurin@donntu.edu.ua, antonina.andrieieva@donntu.edu.ua
Abstract. This paper enlightens the simulation of rock deformation processes
on multicore processing. For mines rock deformation processes affect gas emis-
sions and its prediction is an important research and modeling task. In this pa-
per, the parallel architectures OpenMP and MPI are used. The results show the
modeling time of OpenMP and MPI are better than the serial programming
when the number of grids is given more than 400 points.
Keywords. Parallel simulation, SIMD, MIMD, deformation process, rock de-
formation.
1 Introduction
With the transition of mining operations to the region of great depths (significantly
lower than the gas weathering zone), external and internal factors related to the coal
layer, which cause the occurrence of GDH are likely to change. In this aspect, the
most careful attention should be paid to the geomechanical processes occurring in the
near-bottom part of the coal layer.
In various studies of the stress state of a rock mass, as a rule, the mass is consid-
ered as an elastic medium. The stress concentration in the edge of the layer, located at
depths of more than about 300 - 400 m, leads to the formation of an inelastic deforma-
tion zone [1, 2, 3, 4, 5]. However, taking into account the physical and mechanical
properties of the coal layer, the geomechanical, chemical and gas-dynamic processes
occurring in this zone, besides plastic deformation, brittle destruction of the rock and
deformation of the genetic return also occur. It will be generally fair to consider this
zone as a zone of inelastic deformations.
In the process of coal mining in the bottom, the boundary of the zone of inelastic
deformation moves deep into the massif and the bottom-hole part of the coal layer is
extracted into the mine.
Over time and under the influence of existing stresses, the equilibrium state of the
bottom-hole part of the coal layer can be disturbed, and often it happens. The interac-
tion of these factors, and the physicomechanical properties of the layer and host rocks,
leads to a critical state of coal in the near-bottom part of the layer and its natural sof-
tening and destruction.
It should be noted that the destruction of the edge of the coal layer occurs during
unloading due to excavation of the next strip (part) of coal in the lava. The condition
Copyright © 2020 for this paper by its authors. Use permitted under Creative
Commons License Attribution 4.0 International (CC BY 4.0).
on the destruction of a solids not with an increase or prolonged exposure to normal or
shear stresses, but with unloading, was noted back in the 40-50s of the last century by
the American scientist P. Bridgman [6, 7].
Over time, numerous experimental studies have proved the fundamental impossi-
bility of the origin and occurrence of GDH in the unloaded part of the outburst haz-
ardous layer below the gas weathering zone, at great depths.
The formation and development of the zone of inelastic deformations in the bound-
ary part of the layer always occurs during stress redistribution. This part of the layer is
characterized by reduced tension and, consequently, an increased tendency to gas
recovery - gas evolution. With an increase in the depth of mining operations, and
consequently, an increase in rock pressure, the processes of destruction of the edge
part are undoubtedly intensified.
The conditions on the non-outburst hazard of the edge of the outburst layer formed
the basis of many methods for predicting, preventing GDH and monitoring their ef-
fectiveness, which were included in the new standard of the «Ministry of Industry of
Ukraine» [8].
Substantial growth below the gas weathering zones of the extraction strains was
determined, which allows, in this regard, to produce in it a more productive, con-
trolled extraction of coal from the outburst hazardous layer without causing sudden
outbursts.
On this basis, rock deformation processes affect gas emissions and its prediction is
an important research and modeling task.
2 Mathematical description
To simulate the deformation process, a system of equations of mechanics including
the equations:
• continuity
ui , j 0 (1)
• movement
ji ,i F j u j . (2)
Here:
– the density of the material;
u j – the components of the velocity vector;
j,i – the components of the Cauchy stress tensor;
F j – mass forces.
Here is the density of the material; u j are the components of the velocity vec-
tor; j,i are the components of the Cauchy stress tensor; F j - mass forces; the dot
above means the tiem derivative, the index after the decimal point means the corre-
sponding coordinate derivative. The system of equations is closed by defining rela-
tions that specify the behavior of the medium during deformation. When writing the
defining relations, the decomposition of the stress tensor into the spherical and devia-
tor parts is used:
j ,i j ,i S j ,i (3)
Where:
– average pressure;
S j ,i – stress tensor deviator components:
j,i – Kronecker symbol.
Before the onset of a plastic state, the rates of stress and strain changes are related
by a hypoelastic law:
DS j ,i 1
2 0 ( j ,i nn j ,i ), (4)
Dt 3
DS j ,i
S j ,i S ni nj S nj n, j , (5)
Dt
V
K . (6)
V
where K и – compression and shear modules, respectively. The components
of the strain rate tensor j,i rotational speed tensor components j,i are determined
from the relations:
1
j ,i (u j ,i ui , j ), (7)
2
1
j ,i (u j ,i ui , j ). (8)
2
It is accepted that the deformation consists of elastic e and plastic p parts. For
strain rate, this decomposition is written as:
e p . (9)
Plastic deformation is determined in accordance with the equation of the limiting
surface and plastic potential:
f ( j ,i , jp,i ) 0, (10)
g ( j ,i , jp,i ) 0, (11)
g
d jp,i d , (12)
j ,i
Where:
f – load surface (function) equation;
g – plastic potential;
d - determined during deformation from the condition of plasticity;
j,pi – components of plastic (inelastic) deformation.
Limiting surface (fig. 1) in the area of shear deformation in the pressure range
j 0 described by the equation
f1 ( , ) c (13)
and with pressures 0 – by the equation
( 0 ) 2 2
f 2 ( , ) 2 1 0 . (14)
a2 b
Here (e , ) , c c (e , ) – coefficients of internal friction and adhe-
n n
1
sion, ( S j ,i S i , j / 2) – second stress deviator invariant, t – tear strength, 0 –
2
threshold pressure at which material compaction begins, a 1 0 ,
b c 0 .
The normal to the plastic potential g determines the direction of plas-
tic strain increments, where – is the dilatancy coefficient.
It is accepted that during the development of inelastic deformation, a change in
the surface of the limiting state occurs, as well as damage accumulation, which
macroscopically manifests itself in a change in volume, i.e., in medium dilatancy. The
hardening of the medium is described by the relation
c( p ) c0 [1 h( A( p ) D ( p ))] (15)
1
where h – hardening coefficient, d 2((de j ,i ) 2p
p
/ 2) shear plastic strain
2
rate,
*
– critical deformation, after which degradation of the material predominates.
To account for hardening, a linear relationship A ( ) 2 2 / and quadratic - to
p
take into account softening (accumulation of damage)
D( p ) ( p / ) (16)
The effect of pressure on the limiting strain that the material can withstand before
softening begins is taken into account by the expression:
* 0* (1 w / * ) (17)
where 0 – plastic deformation of the onset of fracture in the absence of compres-
*
sive (restraining) pressure; w and – parameters. This method made it possible to
*
describe the transition from brittle to “viscous” behavior with increasing pressure [9].
Two-dimensional deformation processes are considered for conditions of plane de-
formation when u z 0, xz yz zz 0. .
To simulate the process of rock deformation, the approach [10] was used, which is
based on solving the equations of dynamics of an elastic-plastic medium using an
explicit numerical scheme [11].
Fig. 1. Surface view of the limiting state of rocks
The main difficulty in describing the deformation beyond the elastic limit lies in
the fact that the limiting surface is not fixed, it changes during deformation. On the
loading diagrams, a hardening section is observed, where the increment of plastic
deformation occurs with increasing stresses, followed by stress relief, and destruction
occurs. In the course of deformation, not only the limiting surface changes, but also
the ratio between the increments of the shear and volume parts of plastic deformation,
i.e., the direction of the plastic deformation vector. Thus, parameters describing the
behavior of the rock beyond the elastic limit become accumulated plastic deformation
and pressure functions.
It should be noted that the hardening section on the loading diagram in such media
is very difficult to interpret. Formally, strength is described by two parameters: adhe-
sion and internal friction. However, geometric causes contribute to effective harden-
ing. Dilatancy hardening [12, 13] is especially evident in cramped conditions of de-
formation, when an increase in volume leads to an increase in pressure. This phe-
nomenon is of particular importance in geomechanical processes. It was shown in
[14] that an increase in the effective strength of the sample can be observed even with
a slight decrease in the adhesion.
3 Parallel architecture
There are two different approach in multicore parallel computing - Shared memory
parallelization and Distributed memory parallelization:
Shared memory parallelization(SIMD) - architecture where each processors using
the main memory for communicating of data. OpenMP is a main platform shared
memory architecture. OpenMP is a simple instruction to divide a process into some
subprocess which are called threads. OpenMP using in many programming lan-
guage such as C/C++ and Fortran.
Distributed memory parallelization(MIMD) - architecture where each processors
has own local memory which are called nodes. All nodes connected in one com-
munication network. Each node will execute the code and working with data using
own local memory. Then, final result should be distributed to the main(master)
node to gather all results. To work with this architecture is used MPI (Message
Passing Interface). MPI have implementation in programming languages C/C++
and Fortran and other.
In addition to the two main platforms OpenMP and MPI, there is a hybrid platform -
OpenMPI, which uses the capabilities of both plotforms. In addition to OpenMPI,
there is also a full-fledged language for parallel programming of hybrid architectures -
Unified Parallel C.
Also in OpenMP, in addition to the built-in functions for parallelizing calculations,
the built-in capabilities of the processor for parallelization, the so-called Math Kernel
Library(MKL), are used. MCL can compute operations of addition, subtraction and
multiplication over vectors and matrices in parallel.
The parallel speedup can be obtained by:
T1
S ( p) (18)
Tp
Where:
─ S ( p ) is speedup by using parallel programming using p processor,
─ T1 the CPU time by using single processor,
─ T p CPU time by using p processor.
The parallel efficiency can be obtained by:
S ( p)
E 100% (19)
p
In the process of coal mining in the bottom, the boundary of the zone of inelastic de-
formations moves deep into the massif, which increases the number of calculations
during the simulation, and the acceleration of simulation of this process will make it
possible to obtain data much faster and more accurately due to an increase in the size
grid.
The solution of the system of equations is carried out on the computational grid,
which covers the studied region of the medium or the entire sample, the discretization
of the computational region is carried out. On the constructed computational grid, all
functions are approximated and a system of equations is solved. Of course, difference
schemes that have been successfully used for many years, either for problems in the
deformation of the geomedium, are described in detail in [15, 16]. The standard fea-
tures of numerical calculation are used, which are implemented in original software
packages.
The OpenMP program is simple, for instance in C++, the instruction only add syn-
tax ”#pragma omp parallel for” above the looping codes, then, the single core will
automatically divide the loop into several block.
Fig. 2. Block process of OpenMP architecture
The MPI program in Fig. 3 has two domains of computation. The difference part
from block process of OpenMP is located on looping time. Other processors will be
executed time loop similar to the master.
In MPI, each processors need communication to each other for transferring data
from master to slave. MPI have commands “MPI_Send()” and “MPI_Recv()” to
send and receive data. But MPI (MIMD architecture) needs more communication
time than the SIMD architecture and MPI could not guarantee produces best execu-
tion time.
Fig. 3. Block process of MPI architecture
4 Experiments and results
To obtain simulation results and parallel performance, the following computer
specifications are used.
Table 1. The Computer Specifications
Name Type
Operating system CentOS 6
Processors Intel Xeon X5560 2.8 GHz
RAM 16GB
CPU(s) 16
The performance of MPI and OpenMP will be elaborated. The runtime for MPI
and OpenMP for simulation is done for 10 seconds of the process.
Table 2. Computing time
Number Serial (s) Time (s) of OpenMP Time (s) of MPI
of grid
8 cores 16 cores 8 cores 16 cores
200 27.8096 49.5230 61.9437 36.8255 54.3485
400 94.8357 110.2547 99.3012 82.2647 135.8475
800 370.5342 245.9861 240.1843 237.1254 330.1293
1600 1490.5422 690.5362 580.1232 740.9543 997.3753
Table 2 shows the comparison of modeling time in serial parallel using 8 and 16
cores for each type. Form the Table 2, the modeling time in serial for number of grids
200 and 400 is smaller than the modeling time in parallel except MPI type for 8 cores.
This means the parallel processing is not less efficient since the number of grid points
is small, in this case using single core is already enough. However, using the large
number of points (from 800 and greater), the modeling time in parallel is better than
the serial.
From Table 2, the speedup (18) and efficiency (19) performances of OpenMP and
MPI using 8 and 16 cores can be achieved. The speedup and efficiency profiles are
shown in Fig. 2 and 3 respectively. In Fig. 2, the speedup of OpenMP is better than
MPI for different number of cores. However, the speedup using MPI shown slightly
better in number of grids 800 8 cores. The MPI fails to get the good speedup perform-
ance due to the time for changing messages in processors is dominant than the com-
puting time.
Fig. 4. Speedup of OpenMP and MPI using 8 (left) and 16 (right) cores.
Fig. 5. Efficiency of OpenMP and MPI using 8 (left) and 16 (right) cores.
The efficiency profiles in Fig. 3 are given. The efficiency of parallel processing us-
ing 8 cores is better than 16 cores in both platforms.
5 Conclusion
Two parallel architectures OpenMP and MPI have been implemented for simulat-
ing the 2D deformation process. The results show the modeling time of OpenMP and
MPI are better than the serial programming when the number of grids is more than
200 points.
References
1. Bobrov, I.V., Krichevskiy, R.M., Mihaylov, B.I.: Vnezapnyie vyibrosyi uglya i gaza na
shahtah Donbassa, Ugletehizdat (1954)
2. Nikolin, V.I., Mordasov, V.I., Savchenko, P.I.: Obosnovanie minimalnogo chisla ne-
schastnyih sluchaev, neobhodimyih dlya dostovernogo analiza travmatizma. In: Izvestiya
Gornogo institute, vol. 2, pp. 41-44 (1999)
3. Ruppeneyt, K.V.: Nekotoryie voprosyi mehaniki gornyih porod, Ugletehizdat (1954)
4. Bridgmen, P. W.J.: Appl. Phys (1960)
5. Buzer, G.D., Hiller, K.H., Serdengekti, S.: Vliyanie porovoy zhidkosti na deformatsionnoe
povedenie gornyih porod pri trehosnom szhatii. In: Mehanika gornyih porod, pp. 372-406
(1966)
6. Zabigaylo, V.E., Nikolin, V.I.: Vliyanie katageneza gornyih porod i metamorfizma ugley
na ih vyibrosoopasnost, Naukova dumka (1990)
7. Bridzhmen, P.V.: Issledovanie bolshih plasticheskih deformatsiy razryivov, Izd-vo
inostrannoy literaturyi (1955)
8. Levkin, N.B.: Predotvraschenie avariy i travmatizma v ugolnyih shahtah Ukrainyi, Don-
bass (2002)
9. Stefanov, Yu.P.: Chislennoe modelirovanie deformirovaniya i razrusheniya gornyih porod
na primere rascheta povedeniya obraztsov peschanika. In: Fiziko-tehnicheskie problemyi
razrabotki poleznyih iskopaemyih, vol. 1, pp. 73-83 (2008)
10. Stefanov, Yu.P.: Nekotoryie osobennosti chislennogo modelirovaniya povedeniya uprugo-
hrupkoplastichnyih materialov. In: Fizicheskaya mezomehanika, vol. 8, pp. 129-142
(2005)
11. Uilkins, M.L.: Raschet uprugoplasticheskih techeniy. In: Vyichislitelnyie metodyi v
gidrodinamike, pp. 212-263 (1967)
12. Nikolaevskiy, V.N.: Mehanicheskie svoystva gruntov i teoriya plastichnosti. In: Mehanika
tverdyih deformiruemyih tel. Itogi nauki i tehniki, vol. 6, pp. 5-85 (1972)
13. Garagash, I.A., Nikolaevskiy, V.N.: Neassotsiirovannyie zakonyi techeniya i lokalizatsii
plasticheskoy deformatsii. In: Uspehi mehaniki, vol. 12, pp. 131-183 (1989)
14. Stefanov, Yu.P.: Numerical investigation of deformation localization aтd crack formation
in elastic brittle-plastic materials. In: J.Fract., vol. 128, pp. 345-352 (2004)
15. Wilkins, M.L.: Computer Simulation of Fracture. Lawrence Livermore Laboratory, Rept.
UCRL-75246 (1972)
16. Wilkins, M.L.: Computer Simulation of Dynamic Phenomena, Berlin–Heidelberg–New
York: Springer-Verlag (1999)