Acta Polytechnica CTU Proceedings https://doi.org/10.14311/APP.2023.40.0061 Acta Polytechnica CTU Proceedings 40:61–68, 2023 © 2023 The Author(s). Licensed under a CC-BY 4.0 licence Published by the Czech Technical University in Prague ANALYSIS OF TUNNEL EXCAVATION BASED ON LINEAR DFN-FEM MODELLING Martin Lebeda∗, Petr Kabele Czech Technical University in Prague, Faculty of Civil Engineering, Department of Mechanics, Thákurova 7, 166 29 Prague 6, Czech Republic ∗ corresponding author: martin.lebeda@fsv.cvut.cz Abstract. Simulations of tunnel excavations have to take into account the natural occurrence of joints and faults in the surrounding rock mass, which dominantly control its mechanical response. In this paper, we present work in progress toward 3D finite element analysis of excavation using equivalent rock-mass properties derived from stochastically generated discrete fracture networks (DFNs). The equivalent stiffness is determined by volume averaging. Presently, we solve the problem linearly for an incremental change of the stress state. The fracture’s stiffness is assumed to depend on the initial normal stress acting in direction normal to it. However, within the solved incremental step, we assume the fracture’s stiffness to be constant. This assumption is acceptable for small stress changes. Since the fractures represented in the DFN model have preferred directions, the equivalent stiffness is anisotropic. Keywords: Fractured rock mass, discrete fracture network (DFN), finite element method (FEM), averaging procedure. 1. Introduction According to the current state of art in the field of rock mechanics, brittle structures (or fractures, in- cluding faults and joints) have the dominant effect on the overall mechanical response of the rock mass [1]. Underground structures, such as tunnels, are often situated in rock masses with naturally existing frac- tures. Since it is necessary to take into account the fractures’ influence on the deformation due to the changes in the stress field during excavation, estimat- ing mechanical properties of the rock mass is funda- mental for efficient and safe design of tunnel struc- tures. The finite element method (FEM) is a tool commonly used for simulating tunnel excavation. In the case of a highly fractured rock mass, it would not be effective and feasible to explicitly represent each distinct fracture in the finite element mesh. It is admissible to include only a few largest fractures or localized deformation zones using e.g. interface elements, while the rest of the fracture network is represented by so-called equivalent continuum. The overall effective properties can be determined by aver- aging or homogenization procedure [2]. In this paper, we use the averaging procedure proposed by Oda et al. [3] to determine the effective stress-strain re- lation. As the overall properties of fractured rock mass depend on the size of the simulated area, it is suitable to include the concept of representative volume element (RVE) or statistical volume element (SVE) [4]. The discrete fracture network (DFN) method is one of possible ways to describe the fractures’ geom- etry while taking into account data obtained from structural-geological survey, e.g. [5–8]. One of the most used approaches are the stochastically generated DFNs, which extrapolate the in-situ data by means of statistical probability distributions. In this paper, we present work in progress toward developing and implementing in computer code an ap- proach that connects averaging procedure performed on 3D stochastic DFN with finite element analysis of tunnel excavation. At the present stage of develop- ment, we adopt some simplifying assumptions: (A) We assume that the rock mass in the whole do- main of the finite element (FE) model has uniform equivalent effective mechanical properties. (B) The effective properties are derived by averag- ing [3] over a representative volume element (RVE) of the rock mass, in which shapes, sizes, and ori- entations of individual fractures are modeled by means of stochastic DFN. The RVE size is based on a previous study [9] and corresponds with the size of the FE model. (C) Individual fractures are treated by means of the “parallel plate model” [3], which means that the deformation response of a fracture to normal and shear stress is represented by two parameters ac- counting for the fracture’s roughness, friction angle etc. as well as its size and normal stress acting on it. (D) The FE analysis is performed linearly for an in- cremental change of the stress state. The fracture’s stiffness is assumed to depend on the initial normal stress acting in direction normal to it. However, within the solved incremental step, we assume the stiffness to be constant. We are aware that some of the assumptions may be limiting the validity of the procedure. For example, 61 https://doi.org/10.14311/APP.2023.40.0061 https://creativecommons.org/licenses/by/4.0/ https://www.cvut.cz/en Martin Lebeda, Petr Kabele Acta Polytechnica CTU Proceedings Parameter Value Volumetric density P30 [1/m3] 2.0 Power law exponent α [-] 3.4 Power law minimum fracture size amin [m] 0.3 Fisher distribution concentration factor κ [-] 106 Fisher distribution mean vector [-] variable Volume of the DFN model [m3] 100 × 100 × 100 RVE volume [m3] 15 × 15 × 15 Table 1. Parameters of DFNs. the size of RVE (B) should be chosen not only based on the scatter of the overall properties [9], but also with respect to the gradient of the stress field in the FE model. Or, the linear solution (D) is acceptable only for small stress changes. These limitations will be resolved in a follow-up work. 2. DFN model – description of fractures’ geometry Quantitative information about fracturing of rock mass is often acquired by structural-geological map- ping of fracture traces on “observation windows”, such as rock outcrops or tunnel walls. Typically, the mini- mum recorded trace length is on the order of 10−1 m. Due to the enormous number of corresponding frac- tures in the rock volume, on the one hand, and the limited information available on the observation win- dows, on the other hand, the geometry of the frac- ture network is described by stochastically generated DFN. In this approach, the size and orientation of the fractures are assumed to follow certain probabilis- tic distributions, whose parameters are identified so as to match (in a statistical sense) the traces’ sizes, directions and density within the observation win- dows. A DFN is then generated as a set of spatially distributed polygons, whose sizes and orientations respect the calibrated probabilistic distributions [8]. Since each fracture is explicitly described in the DFN model, the network cannot be efficiently discretized by a 3D FE mesh. However, it provides a suitable basis for the use of an averaging procedure. The fracture network models used in this study were generated by the DFraM software [8]. Fractures’ centers are positioned in the model volume by us- ing Poisson random generation process. Fractures’ sizes are controlled by power law distribution with parameters amin, which is the minimum fracture size (location parameter), and α, which is the law’s ex- ponent (shape parameter). Size of the fracture a is defined as the radius of the circle circumscribed to the fracture, while the fractures are modelled as squares. The orientation of fractures is determined by Fisher distribution, with parameters µ, which is the mean unit normal vector of fractures, and κ, which is the concentration parameter. As the purpose of this study is to demonstrate the solution methodology, we adopt some simplifica- tions while modelling the fracture networks. We use DFNs consisting of one geological set of fractures. Val- ues of the power law distribution parameters are the same for all models and they are based on a previous study [9]. We consider that all fractures are nearly parallel, which is achieved by setting the concentra- tion parameter of Fisher distribution κ to a high value of 106. We will, however, consider fracture sets with different orientations, for which different mean nor- mal vectors µ will be used (see Subsection 6.1). The volumetric fractures density P30, which controls the total number of fractures in the model, was approxi- mately set based on report [10]. Parameters used in this study are listed in Table 1. 3. Overall stiffness of the fractured rock mass evaluated by averaging procedure The calculation of overall effective parameters of frac- tured rock mass is based on the volume averaging relation proposed by Oda et al. [3]: εij = 1 E [ (1 + ν)δikδjl − νδij δkl + ( 1 kn − 1 ks ) Fijkl + 1 4ks (δikFjl + δjkFil + δilFjk + δjlFik) ] σkl =Cijklσkl, (1) where E and ν are Young’s modulus and Poisson’s ratio of the intact rock, respectively, kn and ks are nondimensional parameters related to the fracture’s normal and tangent stiffness, respectively, δij is Kro- necker’s delta, and Fij and Fijkl are so-called second and fourth rank crack tensors, respectively: Fij = 1 V M∑ p=1 S(p)L(p)n (p) i n (p) j , (2) Fijkl = 1 V M∑ p=1 S(p)L(p)n (p) i n (p) j n (p) k n (p) l . (3) Here V is the volume of the domain over which averaging is performed, S(p) is the area of p-th fracture 62 vol. 40/2023 Analysis of tunnel excavation based on linear DFN-FEM modelling Figure 1. Normal-stress vs. fracture-closure dia- gram [12] and schematic indication of the tangent stiffness evaluation. inside the sampling volume, L(p) is a typical size of the fracture equal to a diameter of a circle with the same area and n(p)i are the components of a unit vector normal to the fracture. Note that the typical size L is assumed as a property of the fracture independent on the size of the sampling volume and it is calculated using the original area of the fracture generated in stochastic DFN. The calculation procedure of parameters kn and ks is described in detail in our previous work [11]. In brief, the parameters are expressed as: kn = 3π 8 + κn · L E , (4) ks = 3π 8 + κs · L E , (5) where κn and κs are the normal and shear stiffness (stress/relative displacement) of a fracture and L is the fracture’s size. The normal stiffness parameter κn is to be obtained experimentally. It depends on whether the fracture undergoes loading or unloading and it also varies with pressure acting on it. For the sake of simplicity we determined κn as tangent stiffness by piece-wise linear approximation of the unloading path of the exper- imental normal-stress vs. fracture-closure diagram (Figure 1), which was presented in [12]. The inter- vals of normal stress σn and corresponding tangent stiffness parameters κn and kn are listed in Table 2. Compressive stress is in this study considered with a negative sign. The shear stiffness of fractures is determined using the formula presented in [13] with parameters listed in Table 3: κs = 100 L · |σn| · tan ( J RC · log10 ( J CS |σn| ) + ϕr ) , (6) where J RC is joint roughness coefficient, J CS is joint wall compressive strength and ϕr is residual friction angle. Normal stress |σn| κn kn [MPa] [ MPam ] [-] 0.0–2.5 75 472 1.308 2.5–5.0 210 526 1.540 5.0–7.5 277 778 1.656 7.5–10.0 625 000 2.253 10.0–12.5 781 250 2.522 12.5–15.0 892 857 2.714 15.0–20.0 1 250 000 3.328 20.0–25.0 1 250 000 3.328 25.0–30.0 1 315 789 3.442 30.0–33.0 2 790 698 5.979 Table 2. Fracture’s tangent stiffness (in direction normal to the fracture). Parameter Value Joint roughness coefficient J RC [-] 6.55 Joint wall compressive strength J CS [MPa] 90.0 Residual friction angle ϕr [deg] 20.0 Table 3. Fracture’s shear stiffness parameters. The described constitutive model has been imple- mented by means of Python API as user-defined ma- terial in the open-source FE code OOFEM [14]. It is noted that Equation (3) represents a compliance rela- tion, while in FEM, the stiffness tensor is necessary. Thus, the compliance tensor is first calculated using Equation (3) and then it is numerically inverted to obtain stiffness. 4. Representative volume element As we discussed in the previous section, the overall mechanical properties of the jointed rock mass are evaluated by volume averaging of the fracture net- work. The size of the DFN sub-domain, over which the averaging is performed (V in Equations (2) and (3)), should be chosen large enough to guarantee that the evaluated overall properties are the same for any stochastic realization of the DFN. Then, the sub- domain can be called a representative volume element. Based on a previous study [9], we use RVEs with di- mensions of 15 × 15 × 15 m. For this size the coefficient of variation (COV) of apparent moduli evaluated for 10 realizations of DFN (with similar parameters as here) was less than 20 %. Furthermore, it was shown in [9], that further increasing the sampling volume size did not lead to any significant reduction of the moduli COV. To eliminate potential border effects, the RVEs are obtained by cropping much larger DFN, which is generated in the volume of 100 × 100×100 m. 63 Martin Lebeda, Petr Kabele Acta Polytechnica CTU Proceedings (a). State prior to excavation. (b). State after tunnel excavation. Figure 2. Finite element model. 5. 3D FE modelling of tunnel excavation To demonstrate the use of the DFN-based constitutive model, we numerically simulate by FEM the excava- tion of a fictitious tunnel in rock at the depth of 550 m. The tunnel has a circular cross-section with the di- ameter of 3.0 m. The tunnel is modelled in a cubic body, representing the surrounding rock, with side size of 15 m (Figure 2). Thus, the distance between the tunnel and the boundary of the model is at least twice the diameter of the construction. The coordinate system is introduced with the ver- tical axis z oriented upwards. The edges of model are parallel with the coordinate axes and the tunnel is parallel with the y-axis. The finite element mesh consists of 3D hexahedrons and it was created with Salome software [15]. The simulation is performed in two steps. The first step (time t = 0) represents the initial state prior to excavation – Figure 2a. The displacement of the rock body is constrained by statically determinate supports and surface tractions, corresponding to the initial geostatic pressure, are applied on all sides of the body. The initial stress state is idealized as uniform within the modeled domain. The initial stress tensor is based on field measurements from Rožná mine [16]. Using data from vertical boreholes it was found that the directions of the principal stresses are vertical and horizontal. The vertical principal stress is calculated by Equation (7): Sv = −ρ · g · h, (7) where ρ is volumetric weight, g is gravity accelera- tion and h is depth beneath the ground level. The parameters for calculation of the vertical principal stress (Table 4) are adopted from [17]. As the major and minor horizontal principal stresses, SH and Sh, respectively, we use the average of the values mea- sured on different levels of the mine. The values of the Parameter Value Major horizontal principal stress SH [MPa] −22.8 Minor horizontal principal stress Sh [MPa] −14.6 Volumetric weight ρ [ kg 3 m ] 28 000 Depth of the mine beneath the ground level [m] 550 Vertical principal stress Sv [MPa] −15.4 Table 4. The initial stress state. initial stress are listed in Table 4. Orientation of the model is set so, that the tunnel axis is parallel with the minor horizontal principal stress, which means that the direction of Sh corresponds to the global axis y. It should be noted that, even though the stress field is uniform throughout the analyzed domain in the first step, the normal traction σn acting on each fracture in the DFN is generally different due its different direction. This traction, in turn, affects the fractures’ normal and shear stiffness, as seen in Equation (6) and Table 2, which must be taken into account in the averaging process. In the second calculational step (time t = 1), the excavation of the tunnel is simulated by removing the finite elements inside the tunnel space while keeping the boundary conditions on the outer surface of the rock body unchanged. This creates a new traction- free boundary on the tunnel walls, which results in deformation and change of stress state in the rock. We should remark that, at the present stage of progress, the problem is solved as incrementally linear. That is, the second step is solved with constant tangent stiffnesses of the fractures, which correspond to the initial stress state. 64 vol. 40/2023 Analysis of tunnel excavation based on linear DFN-FEM modelling Figure 3. Horizontal convergence ∆ux calculated for different DFN models. Figure 4. Vertical convergence ∆uz calculated for different DFN models. 6. Case examples 6.1. Fracture network models To demonstrate the proposed simulation procedure of 3D underground structure excavation, we use sev- eral simple examples. In all of the examples, the fracture network consists of a set of nearly parallel fractures with the same DFN parameters, except the mean normal vector. In particular, the simulations are performed on four types of DFN models, in which fractures are: perpendicular to axis x/y/z and at an angle of 45 degrees between negative x axis and positive z axis, see Table 1. Even though the size of the RVE for averaging of the rock properties was rationally selected, 3 different stochastic realizations of DFN models with different random seed value were generated for verification. Thus, in total we work with 12 DFN models. 6.2. Results and discussion 6.2.1. Displacements Figure 3 and Figure 4 present the calculated values of the horizontal and vertical convergence for each DFN case. The tunnel convergence is evaluated as the relative displacement between points on the tunnel walls caused by the tunnel excavation (that is, the incremental displacement between the final displace- ment at time t = 1 and the initial displacement at t = 0). Positive values of the convergence correspond to the closing of the tunnel walls. The points at which Figure 5. Location of reference points and lines used to report results. the displacement was recorded are pictured in Fig- ure 5. Figure 6 through Figure 9 show the incremental displacement of the model boundaries including the excavated tunnel. Considering that the initial vertical pressure Sv is about 2⁄3 of the major horizontal pressure SH and 65 Martin Lebeda, Petr Kabele Acta Polytechnica CTU Proceedings that SH acts perpendicularly to the tunnel axis, the displacements of the tunnel walls qualitatively corre- spond with expectations: • For a fracture set perpendicular to tunnel axis, i.e. perpendicular to the axis y, the effect of the frac- tures on the convergence is rather small and the convergence is the least among all the analyzed cases (the second group of bars in Figure 3 and Fig- ure 4). Consistently with the initial stress state, the vertical convergence is smaller than the horizontal one, which is also obvious from the deformed shape of the tunnel cross-section (Figure 6). • Referring to the cases with vertical and horizontal fractures parallel with the tunnel axis (i.e. per- pendicular to axis x and perpendicular to axis z, respectively), it is seen in Figure 3 and Figure 4 that the convergence in the direction normal to the fractures significantly increases in comparison with the previous case. The convergence as well as the shape of the cross-sectional distortion reflect the fractures-induced anisotropy. In the rock mass with vertical fractures, the higher initial horizontal stress results in horizontal convergence, which is signifi- cantly larger than the vertical one (Figure 7). On the other hand, when the fractures are horizontal, the resulting horizontal and vertical convergences are close to each other (Figure 8). • Finally, when the fractures are inclined at 45 de- grees, the horizontal and vertical convergences fall between those observed with the vertical and hor- izontal sets. The axes of the distorted shape of the tunnel (Figure 9) are close to, but do not align with the fractures, which can be attributed to the difference of the initial horizonal and vertical stress and normal and tangential displacement occurring on the fractures. Figure 3 and Figure 4 also indicate that the larger components of convergence, which are strongly af- fected by the fractures, exhibit notable variation among the different stochastic realizations of the DFN. This finding suggests that the employed sampling vol- ume was still not large enough to meet the criteria as an RVE. 6.2.2. Stress Figure 10 and Figure 11 show the distribution of nor- mal stresses along the horizontal and vertical lines indicated in Figure 5. The results are reported for one realization of the DFN model with fractures at the an- gle of 45 degrees. Both plots show the uniform initial stress state at time t = 0 and the stress distribution after excavation at t = 1. The difference between the lines corresponds to the incremental change of stress due to excavation. The plots show that, as the tunnel is excavated and a new traction-free boundary is introduced at the tunnel wall, the radial pressure correctly decreases from the initial uniform state and tends to zero at the tunnel wall. On the other hand, (a). Vertical section through the DFN. (b). Deformed shape of the tunnel (red). Figure 6. Model with fractures perpendicular to y axis. (a). Vertical section through the DFN. (b). Deformed shape of the tunnel (red). Figure 7. Model with fractures perpendicular to x axis. (a). Vertical section through the DFN. (b). Deformed shape of the tunnel (red). Figure 8. Model with fractures perpendicular to z axis. (a). Vertical section through the DFN. (b). Deformed shape of the tunnel (red). Figure 9. Model with fractures at the angle of 45 de- grees between negative x axis and positive z axis. 66 vol. 40/2023 Analysis of tunnel excavation based on linear DFN-FEM modelling Figure 10. Model with fractures at the angle of 45 degrees: distribution of normal stresses σx and σz along horizontal line A-A’ (Figure 5). the pressure in the circumferential direction increases, which is also to be expected. Note that fractures are distributed at an angle around the tunnel, so some experience an increase of pressure (loading) while pres- sure on others decreases (they unload). The stress increment on fractures in the close vicinity of the tunnel may be also quite large. Therefore, for some fractures, using the constant tangent stiffness esti- mated from the unloading path of the stress-closure diagram (see Section 3 and Figure 1) might have not been accurate. 7. Conclusions The 3D FEM simulation of tunnel excavation with rock mass parameters evaluated by averaging proce- dure based on DFN has been examined. This paper has been focused on the demonstration of the solu- tion methodology. The results allow us to draw the following conclusions: • Deformed shapes calculated by FE simulations are consistent with the orientation of fractures in DFNs and the initial stress state. Largest tunnel conver- gence is observed in the directions close to the least overall stiffnesses of the fractured rock mass. • Some difference in the displacements calculated with different stochastic realizations of the DFNs was observed. The appropriate size of the sampling volume should be further investigated and if RVE cannot be found, the concept of SVE might be adopted. • The largest changes in the stress state after the tunnel excavation were observed close to tunnel, Figure 11. Model with fractures at the angle of 45 degrees: distribution of normal stresses σx and σz along vertical line B-B’ (Figure 5). which is in agreement with the expectations. On the other hand, on the outer boundaries of the models, the stresses after excavation were slightly different from the initial state, which indicates that a larger domain of the model is required to eliminate the boundary effect. • Considering the observed stress variations after the tunnel excavation, an incremental nonlinear calcu- lation with fractures’ stiffness updated according to the actual stress state and loading-unloading may be necessary for a more accurate prediction of the tunnel convergence. • Although it was not elaborated in the present work, the demonstrated analysis method provides infor- mation about the changes of the pressure acting on the fractures, which affect their hydrogeological properties. The proposed methodology, therefore, may provide valuable inputs for hydrogeological simulations of ground water flow and contaminants’ transport. Acknowledgements This paper was financially supported by Czech Technical University in Prague under SGS project no. SGS22/030/OHK1/1T/11. References [1] L. Jing. A review of techniques, advances and outstanding issues in numerical modelling for rock mechanics and rock engineering. International Journal of Rock Mechanics and Mining Sciences 40(3):283–353, 2003. https://doi.org/10.1016/S1365-1609(03)00013-3 67 https://doi.org/10.1016/S1365-1609(03)00013-3 Martin Lebeda, Petr Kabele Acta Polytechnica CTU Proceedings [2] M. Hori, S. Nemat-Nasser. On two micromechanics theories for determining micro–macro relations in heterogeneous solids. Mechanics of Materials 31(10):667–682, 1999. https://doi.org/10.1016/S0167-6636(99)00020-4 [3] M. Oda, T. Yanabe, Y. Ishizuka, et al. Elastic stress and strain in jointed rock masses by means of crack tensor analysis. Rock Mechanics and Rock Engineering 26(2):89– 112, 1993. https://doi.org/10.1007/BF01023618 [4] X. Yin, W. Chen, A. To, et al. Statistical volume element method for predicting microstructure–constitutive property relations. Computer Methods in Applied Mechanics and Engineering 197(43-44):3516–3529, 2008. https://doi.org/10.1016/j.cma.2008.01.008 [5] W. S. Dershowitz, P. Wallmann, S. Kindred. Discrete fracture modelling for the Stripa site characterization and validation drift inflow predictions, STRIPA-TR-91-16, 1991. SKB Stockholm. https://www.researchgate.net/profile/William_ Dershowitz/publication/266334817_Discrete_ Fracture_Modelling_for_the_Stripa_Site_ Characterization_and_Validation_Drift_Inflow_ Predictions/links/563bb3ea08ae34e98c47cfaf.pdf [6] R. Munier. Statistical analysis of fracture data, adapted for modelling Discrete Fracture Networks-Version 2, SKB R-04-66, 2004. SKB - SvenskKärnbränslehantering AB Swedish Nuclear Fuel and Waste Management Co. http: //skb.se/upload/publications/pdf/R-04-66.pdf [7] C. Darcel, P. Davy, O. Bour, J.-R. D. Dreuzyr. Discrete fracture network for the Forsmark site, SKB Rapport R-06-79, 2006. SKB - SvenskKärnbränslehantering AB Swedish Nuclear Fuel and Waste Management Co. http://www.iaea.org/inis/collection/ NCLCollectionStore/_Public/38/013/38013872.pdf [8] P. Kabele, et al. Mathematical Modeling of Brittle Fractures in Rock Mass by Means of the DFN Method – Final report SÚRAO TZ 286/2018/ENG, 2017. SÚRAO - Radioactive Waste Repository Authority. [9] M. Lebeda, P. Kabele. The effect of sampling volume size on the apparent stiffness of jointed rock mass. APP 34:38–42, 2022. https://doi.org/10.14311/APP.2022.34.0038 [10] L. Gvoždík, et al. Transport of radionuclides from deep geological repository/Testing of conceptual and numeric models – Final report SÚRAO TZ 463/2020/ENG, 2020. SÚRAO - Radioactive Waste Repository Authority. [11] M. Lebeda. Určení tuhosti porušeného horninového masivu s využítím DFN modelů. Master’s thesis, Czech Technical University in Prague, Faculty of Civil Engineering, Department of Mechanics, 2020. [12] S. C. Bandis, A. C. Lumsden, N. R. Barton. Fundamentals of rock joint deformation. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts 20(6):249–268, 1983. https://doi.org/10.1016/0148-9062(83)90595-8 [13] N. Barton, V. Choubey. The shear strength of rock joints in theory and practice. Rock Mechanics 10(1-2):1– 54, 1977. https://doi.org/10.1007/BF01261801 [14] B. Patzák. OOFEM project home page, 2000. [2022-10-13]. http://www.oofem.org [15] EDF, CEA. Salome, 2001. [2022-10-17]. https://www.salome-platform.org/ [16] Z. Bukovská, a kol. Získání dat z hlubokých horizontů dolu Rožná, TZ 464/2020, 2020. SÚRAO - Správa úložišť radioaktivních odpadů. [17] K. Souček, M. Vavro, L. Staš, et al. Komplexní geologická charakterizace prostorů PVP Bukov – část II geotechnická charakterizace, Závěrečná zpráva SÚRAO ZZ 221/2018, 2017. SÚRAO - Správa úložišť radioaktivních odpadů. 68 https://doi.org/10.1016/S0167-6636(99)00020-4 https://doi.org/10.1007/BF01023618 https://doi.org/10.1016/j.cma.2008.01.008 https://www.researchgate.net/profile/William_Dershowitz/publication/266334817_Discrete_Fracture_Modelling_for_the_Stripa_Site_Characterization_and_Validation_Drift_Inflow_Predictions/links/563bb3ea08ae34e98c47cfaf.pdf https://www.researchgate.net/profile/William_Dershowitz/publication/266334817_Discrete_Fracture_Modelling_for_the_Stripa_Site_Characterization_and_Validation_Drift_Inflow_Predictions/links/563bb3ea08ae34e98c47cfaf.pdf https://www.researchgate.net/profile/William_Dershowitz/publication/266334817_Discrete_Fracture_Modelling_for_the_Stripa_Site_Characterization_and_Validation_Drift_Inflow_Predictions/links/563bb3ea08ae34e98c47cfaf.pdf https://www.researchgate.net/profile/William_Dershowitz/publication/266334817_Discrete_Fracture_Modelling_for_the_Stripa_Site_Characterization_and_Validation_Drift_Inflow_Predictions/links/563bb3ea08ae34e98c47cfaf.pdf https://www.researchgate.net/profile/William_Dershowitz/publication/266334817_Discrete_Fracture_Modelling_for_the_Stripa_Site_Characterization_and_Validation_Drift_Inflow_Predictions/links/563bb3ea08ae34e98c47cfaf.pdf http://skb.se/upload/publications/pdf/R-04-66.pdf http://skb.se/upload/publications/pdf/R-04-66.pdf http://www.iaea.org/inis/collection/NCLCollectionStore/_Public/38/013/38013872.pdf http://www.iaea.org/inis/collection/NCLCollectionStore/_Public/38/013/38013872.pdf https://doi.org/10.14311/APP.2022.34.0038 https://doi.org/10.1016/0148-9062(83)90595-8 https://doi.org/10.1007/BF01261801 http://www.oofem.org https://www.salome-platform.org/ Acta Polytechnica CTU Proceedings 40:61–68, 2023 1 Introduction 2 DFN model – description of fractures' geometry 3 Overall stiffness of the fractured rock mass evaluated by averaging procedure 4 Representative volume element 5 3D FE modelling of tunnel excavation 6 Case examples 6.1 Fracture network models 6.2 Results and discussion 6.2.1 Displacements 6.2.2 Stress 7 Conclusions Acknowledgements References