LSM-DFN Modeling for Seismic Responses in Complex Fractured Media: Comparison of Static and Dynamic Elastic Moduli
LLSM-DFN Modeling for Seismic Responses in Complex Fractured Media: Comparison of Static and Dynamic Elastic Moduli
Ning Liu气Li-YunFu b ·** a college of Mechanical and Electrical Engineering, Beijing University of Chemical Technology, Beijing 100029, China b Key Laboratory of Deep Oil and Gas, China University of Petroleum (East China), Qingdao 266580, Shandong, China ARTICLE I NFO Keywords: LSM-DFN modeling se1s1111c responses complex fractured media static and dynamic elastic moduli
ABSTRACT Crack microgeometries pose a paramount influence on effective elastic characteristics and sonic responses. Geophysical exploration based on seismic methods are widely used to assess and understand the presence of fractures. Numerical simulation as a promising way for this issue, still faces some challenges. With the rapid development of computers and computational techniques, discrete-based numerical approaches with desirable properties have been increasingly developed, but have not yet extensively applied to seismic response simulation for complex fractured media. For this purpose, we apply the coupled LSM-DFN model (Liu and Fu, 2020b) to examining the validity in emulating elastic wave propagation and scattering in naturally-fractured media. By comparing to the theoretical values, the implement of the schema is validated with input parameters optimization. Moreover, dynamic elastic moduli from seismic responses are calculated and compared with static ones from quasi-static loading of uniaxial compression tests. Numerical results are consistent with the tendency of theoretical predictions and available experimental data. It shows the potential for reproducing the seismic responses in complex fractured media and quantitatively investigating the correlations and differences between static and dynamic elastic moduli. Natural fractured reservoirs have long been an important target for the oil and gas industry. Fracture networks ranging from the microscopic to regional scales provide highly permeable conduits for fluid flow and pose a paramount influence on the mechanical and transport properties of rocks. Therefore, understanding the distribution and development characteristics of fractures at different scales is the primary factor determining the success or failure of fractured reservoir exploration and development. To determine the presence of fractures in the subsurface, geophysical exploration with seismic methods or sonic logs measurements are widely used, based on the sensitivity of wave velocities, amplitudes and spectral characteristics to fracture compliance (Vlastos, Liu, Main and Narteau, 2007). While, due to the difficulty that the distribution of fractures is complex with a wide range of spatial scales, and the density is highly variable in space, the geometrical and statistical models of fracture patterns are still debated. Many theoretical models have been proposed to analytically investigate the velocity anisotropy (including shear wave splitting) and scattering attenuation caused by constituent fractures (e.g., Crampin, 1978, 1984; Hudson, 1981, 1986), but limited to the idealization and oversimplification of complex fractures with long wave approximation for low crack densities (Vlastos et al., 2007). Nowadays, numerical method as an efficient alternative for general applicability are increasingly applied to emulate seismic responses in complex fractured media. Numerical simulations as an efficient supplement to experimental measurements, provide an independent verification of theoretical predictions. Generally, numerical approaches are divided into two categories: Continuum methods and discrete-based models (Suiker, Metrikine and De Borst, 2001; Liu, Fu, Tang, Kong and Xu, 2020b). Conventional continuum mechanics takes matter to be continuously distributed throughout a body. It provides a reasonable assumption for analyzing the macroscale behavior of rocks by homogenization of microstructural effects. The latter regards rocks as an assembly of microstructural elements that interact with each other by microstructural forces. For fractured medium modeling, a numerical scheme with discrete fracture networks (DFNs) is increasingly used instead of the equivalent medium theories. Early continuum-based methods, e.g., the finite element method (FEM), finite *Corresponding author 蒋 Principalcorresponding author由nicolali 呻 buaa.edu. en (N. Liu); lfu©upc. edu. en (L. Fu)ORCID(s): 0000-0001-8692-8405 (L. Fu) Page 1 of 18 SM-DFN modeling for seismic responses difference method (FDM), and boundar y element method (BEM), model the macroscopic behavior of fractured media b y implementing the corresponding macroscopic constitutive models (Swoboda, Shen and Rosas, 1998; Zhao, 2015, e.g.,). Then, Lei, Latham and Tsang (2017) combine continuum-based approaches with discrete fracture network (DFN) to model fractured rocks with only a few or plent y of fractures associated with onl y a small amount of displacement/rotation b y introducing interface elements, or joint elements (Goodman, Ta y lor and Brekke, 1968; Lei et al., 2017). However, modeling the high-densit y and complex DFNs remains difficult, which is regarded as the intrinsic limit of continuum-based methods (Jiang, Zhao and Khalili, 2017). For more complex DFNs, discrete-based approaches seem more suitable (Jing and Hudson, 2002), especiall y for fractured rocks with a wide range of mineral compositions and fabric anisotropies (Liu and Fu, 2020a). Discontinuum-based methods, like molecular d y namics (MD), lattice spring model (LSM) (Hrennikoff, 1941), and discrete element method (DEM) (Cundall, 1971), have been widel y used to emulate mechanical deformations in rocks (Zhao, 2010; Liu and Fu, 2020a). Inhomogeneous effects at the microlevel could be captured b y these discretebased approaches (e.g., Suiker et al., 2001; Liu et al., 2020b), where granular textures, particle-scale kinematics, and force transmission are correlated. Harthong, Scholt七sandDonze (2012) propose a coupled DEM-DFN model for strength characterization of rock masses. Then, Bonilla-Sierra, Scholtes, Donze and Elmouttie (2015) anal y ze rock slope stabilit y using DEM-DFN modeling. Based on this modeling approach, Wang and Cai (2019) emulate excavation responses in jointed rock masses. To the best knowledge of the authors, those works mainl y focus on the mechanical properties and crack opening, propagation, or coalescence under static, or d y namic loadings. However, only a few works appl y the discrete-based methods to emulate the wave propagation, (e.g., Suiker et al., 2001; O'Brien and Bean, 2004; Zhu, 2017; Liu et al., 2020b; Cheng, Luding, Saitoh and Magnanimo, 2020), not to mention the issues related to complex fractured media. The aim of this paper is to attempt to combine a discrete-based numerical method with DFN s for the simulation of the elastic wave propagation and scattering in complex fractured media. Among these discrete-based methods, LSM attracts the most interest because it is flexible to model both continuous and discontinuous s y stems in a discrete way (Ostoja-Starzewski, 2002), and can avoid singularit y -related issues in continuum-based numerical simulation methods (Pan, Ma, Wang and Chen, 2018). Moreover, unlike the DEM where elements interact through contact surfaces, the LSM connects elements b y springs or beams, with several desirable properties, such as broad applicability, eas y implementation, and high flexibilit y to handle the contact complexit y of granular materials (Liu et al., 2020b). Recentl y , more advanced LSMs are developed to avoid the Poisson's ratio limitation of the early LSMs (Hrennikoff, 1941), e.g., Born springs (Hassold and Srolovitz, 1989), multi-bod y shear springs (Monette and Anderson, 1994), noncentral sheartype springs (Griffiths and Mustoe, 2001), beam element model (Karihaloo, Shao and Xiao, 2003; Lilliu and van Mier, 2003), distinct lattice spring model (DLSM) (Zhao, Fang and Zhao, 2011), and modified LSM (Liu et al., 2020b). Based on the modified LSM (Liu et al., 2020b), Liu and Fu (2020b) develop a coupled LSM-DFN model to investigate the stress-orientation effect on the effective elastic anisotrop y of complex fractured media. In this stud y , we introduce the coupled LSM-DFN model (Liu and Fu, 2020b; Liu, 2020) to examining the validity in emulating elastic wave propagation and scattering in naturall y -fractured media. B y comparison to the theoretical predictions, the implement of the modified LSM for seismic responses are validated with input parameters optimization. Base on this numerical approach, the d y namic elastic moduli are calculated and compared with static ones from quasistatic loading of uniaxial compression tests. The remaining parts of this paper are organized as follows: The methods used are introduced briefl y , including the discrete-based numerical approach for seismic responses in complex fractured media, an integrated workflow of LSM-DFN modeling, and the methodolog y for estimating static and d y namic elastic moduli from numerical experiments in Sect. 2; to enable high-accurate modeling, the implementation of the modified LSM for seismic responses are validated with "optimal" input parameters b y calibration process in Sect. 3; Section 4 presents the numerical results of LSM-DFN models with a single crack, uniformly oriented cracks, and natural fracture networks, and comparisons are made between the static and d y namic elastic moduli for the fractured samples; The conclusions and future work are underlined in Sect. 5. 2.Methodology2.1. Numerical method: Integrated LSM-DFN scheme To build an LSM-DFN coupling model for fractured media, homogenous matrix and fracture geometries are represented b y lattice spring models (LSMs) and discrete fracture networks (DFNs), respectively. Lattice spring modeling (LSM) methods can be traced back to the 1940s, proposed b y Hrennikoff (1941) to solve the elasticit y problem. Then, Page 2 of 18 SM-DFN modeling for seismic responses ■■■■■■■■■■■■■■■■ .. . c 'l' 吓 a X 二 e � 乍 c x ....... .. .. ··伽 .... fu ... •• .... .. ..... : ... a y Figure 1: Schematic definition of the elastic constants in two dimensions. the lattice model is introduced to calculate the effective mechanical properties (Nayfeh and Hefzy, 1978), distributed disorder influence (Curtin and Scher, 1990a), stress concentrations and toughness increases (Curtin and Scher, 1990b). With the rapid development of computers and computational techniques, LSM is applied to the study for other mechanical and dynamic behaviors, e.g., tectonic processes (Saltzer and Pollard, 1992; Mora and Place, 1998), and fracture problems (Schnaid, Spinelli, Iturrioz and Rocha, 2004; Kosteski, D'ambra and Iturrioz, 2009). To overcome the fixed Poisson's ratio limitation existing in conventional LSMs, some advanced models are developed. For example, Liu et al. (2020b) modify the LSM by introducing an independent micro-rotational inertia to avoid the Poisson's ratio limitation and improve the LSM in characterizing scale-dependent effects. Then, Liu and Fu (2020a) use this modified LSM to study the elastic characteristics of digital cores, to incorporate with random void model (RVM) for organic-rich shale modeling (Liu, Fu, Cao and Liu, 2020a), develop a coupled LSM-DFN model to investigate the stress-orientation effect on the effective elastic anisotropy of complex fractures (Liu and Fu, 2020b) and to consider rough contact deformation for fracture surfaces (Liu, 2020). While, all the works remain quasi-static elastic issues. Here, we introduce the LSM-DFN modeling to emulate elastic wave propagation and scattering in naturallyfractured media. Firstly, we emulate the matrix as a linearly elastic, homogeneous, and isotropic material using the regular triangular modified LSM. For a digital image of fractured rocks, the lattice nodes are set in an approximate zone represented by digital pixels. DFNs can be stochastically created (Hyman, Gable, Painter and Makedonska, 2014), or extracted from a digital image (Liu and Fu, 2020b; Liu, 2020). The resulting DFNs are incorporated into the LSMs to form the LSM-DFN model for fractured rocks. According to stress-free assumption for crack surfaces (Murai, 2007), the corresponding nodal interactions crossing the DNFs are removed from the model (Lisjak and Grasselli, 2014; Liu and Fu, 2020b). 2.2. Estimation of static elastic modulus: Uniaxial compression test The determination of static elastic modulus is based on the measurement of the relation between elastic deformation and the known force (Ciccotti and Mulargia, 2004). For a two-dimensional (2D) case, three independent constants characterize the elasticity of the transversely isotropic model in principal coordinates, as schematically shown in Figure 1.the constitutive relations can be given by, [::,]�
Vxy百° l-丘 vyx-凡1-乌0 y OOl一如 UJ (1) Here, we conduct uniaxial compression tests to calculate the stress-strain responses of specimens and obtain the Young's moduli and the pair of Poisson's ratios in Eq. (1) by the slopes of the response curves, 凡 =
E X (2) E = 5'.. y E ' X (3) Page 3 of 18 nd E V = _ _2'. xy E' X E V = _2_. yx E' y EV = EV x yx y xy• LSM-DFN modeling for seismic responses (4) (5) (6) owing to the known symmetry condition. It is worth noting that Saint Venant (1863) provides an approximation of the shear moduli in orthotropic materials, based on the values of the rest of the elastic constants, 1 1 1 + j i —=—+ µii Ei E (i,j = x,y). 2.3. Estimation of dynamic modulus: Seismic responses (7) Dynamic elastic moduli can be derived from elastic-wave velocity and rock density. Based on the seismic responses, we can measure velocity of a stress wave passing through a material. For carrying out this seismic pulsing method, a vertical downward displacement excitation, Ricker wavelet, given by Uz = -A 。 (1- 2兀千(t-T) 汀 exp(-兀于(t-T) )' (8) where A 。, T and f represent amplitude, period, and central frequency, respectively, is applied at the center of the upper surface, and the receivers are arranged as shown in Figure 2. The velocities of the compressional and shear propagation of elastic waves covering the certain distance for a known receiver can be calculated by the time between the sending and receiving waves. According to the classical theory of elastic wave propagation, the compressional and the shear wave velocity, cp and cs are given by, Cp =尸, cs = l , where pis the density. The Young's modulus in 2D can be obtained by E = 灶 2µ 'so we could estimate the dynamic Young's modulus by, Edyoami ,�(cs) (2- (气 )} Validation and Calibration (9) (10) (11) (12)
The objective of this section is to validate the implement of the modified LSM for seismic responses. By calibrationprocess, this LSM-DFN model could enable high-accurate modeling by properly choosing reasonable input parameters, e.g., meshing resolutions, stability conditions, and numerical damping (Liu and Fu, 2020b). Poisson's ratioµis a keyPage 4 of 18
SM-DFN modeling for seismic responses 士 l “ri Figure 2: Source-receiver placement of seismic response simulation: The source and the receivers are marked by red star and yellow triangles. Table 1 Material properties in the simulation
Parameters Density p (kg/m ) Value
Elastic constants
E (GPa) µ p (m/s) Vs (m/s) ; K =�and the com ress1onal to shear veloc1t Kn p y 上, b C5 y 从 -K v= = + K +{ , (13) and (干玑 +Ks= 二 Cs 氐+K +( (14) More details are given in our previous works (Liu et al., 2020b). We perform a comprehensive series of numerical simulations to explore the "optimal" inputs for the subsequent studies. The physical and mechanical properties of lattice nodes are reported in the Table 1. 3.1. Lattice spacing and meshing sensitivity Lattice arrangement of discrete-based numerical methods may cause numerical anisotropy. Liu and Fu (2020b) study the anisotropic effect and meshing sensitivity of the LSMs by Brazilian tests. With a higher resolution, the anisotropy caused by lattice orientation could be minimized, but with heavy computational costs. In that work, Liu and Fu (2020b) give a suggestion for meshing resolution of a round specimen with a diameter of D and lattice spacing of d: The maximum errors of models are 1.09% for D/d = 400, and 0.22% for D/d = 600. Therefore, we build models with a resolution of D / d = 600 for more accurate simulations. 3.2. Time step and source frequency The central difference method is conditionally stable. To ensure the stability of this explicit integration algorithm applied in the modified LSM, the critical time step ll.tcr is employed by _.if�_.if CT ' (15) Page 5 of 18
T 2T 4T 6T 8T 10T 12T-1.5-1.0-0.50.00.51.0
Time10 t 50 t20 t 100 t Ricker wavelet curve .00 0.20 0.40 0.60 0.80 1.00-4.000.004.008.0012.0016.0020.00
Normalized time=0 =0.02 =0.04=0.06 =0.08 =0.10=0.20
SM-DFN modeling for seismic responses Figure 5: Snapshots of the vertical displacement wavefield using modified LSM at t = 0.10, 0.20, 0.50, 1.00 s Figure 6: Uniaxial compression tests of LSM-DFN models with six different half-lengths of 5, 20, 50, 100, 150, and 200 m and seven different inclined angles of 0 ° , 15 ° , 30 ° , 45 ° , 60 ° , 75 ° , and 90 ° . 4.1. Case 1: LSM-DFN with a single crack In this case, we couple a single crack located at the center of the domain into the lattice spring networks. The single crack is designed with six different half-lengths of 5, 20, 50, 100, 150, and 200 m and seven different inclined angles of 0 ° , 15 ° , 30 ° , 45 ° , 60 ° , 75 ° , and 90 ° . Totally, 6 x 7 X 2 = 84 numencal expenments are conduct for quasi-static loadings and seismic responses. Figure 6 presents the loading patterns under uniaxial compression for static moduli. We hold the LSM-DFN model with a half-length of 200 m up as an example shown in Figure 7: From left to right, we present the differential stress fields under quasi-static loading state at Ey = 0.60%, and snapshots of the vertical displacement wavefield at t = 0.7 s. From this figure, we can see static stress concentration appears at the tips of these cracks under uniaxial compression tests as theoretically predicted; a large fracture hinders the wave propagation and compressional wave are reflected upwardly at the crack interface for seismic responses. For a two-dimensional homogeneous and isotropic material containing a single crack, the effective Young's modulus E in the direction of the loading is given by Walsh (1965), ½� 土(, + 兀 �'cos p), (24) where Eintact• L e , /3, and A is the effective Young's modulus of the intact sample, the half-length of a crack, the inclined angle with the horizonal axis, and the area of a specimen, respectively. Figure 8 shows the counter plots of the effective Page 8 of 18 .00 0.04 0.08 0.120306090 Crack density (-) 0.000.2000.4000.6000.8001.00E t /E intact I n c li n e d a ng l e ( ° ) Crack density (-) 0.000.2000.4000.6000.8001.00E s /E intact d /E intact
50 100 150 2000.00.20.40.60.81.01.2-1 0 1 2 3 40.960.970.980.991.001.01 intact =1.0 (L c ) x (m)Theoretical solution SM-DFN modeling for seismic responses (a)L e= (b)L e= Differential stress I 劲p GPa) 。 Differential stress ooooo 户颌p~ (c)L e= Figure 10: Differential stress fields under quasi-static loading state at EY = 0.60% of LSM-DFN models with crack halflengths of 50 m, 100 m, 200 m and four different inclined angles of 30 ° , 45 ° , 60 ° , and 90 ° . of the perturbations caused by all the secondary waves. The phase of each secondary wave is delayed by a quantity corresponding to the distance from the center of the receiver, so the arriving secondary waves interfere with each other, resulting in diffraction (Pao and Mow, 1971). Figure 13 shows a comparison made between the theoretical, static, and dynamic moduli. It can be seen from the figure that the change trend of the static elastic moduli is in good agreement with the theoretical values. The models with the first and second crack arrangements can be seen based on the Voigt and Reuss bound theories, respectively. The cracked part can be regarded as a material layer with a lower stiffness. The effective moduli of the models with the third spatial arrangement are within the zone covered by Voigt-Reuss bounds as theoretical predictions. These numerical experimental results further verify the effectiveness of the LSM-DFN models. However, due to stress stiffening caused by uniaxial compression, the static Young's modulus values are slightly higher than the theoretically predicted ones. However, for Case 2, there are obvious deviations between the elastic moduli obtained from the seismic responses and the theoretical predictions, especially for the dynamic elastic modulus calculated from Receive 4. The phenomenon is related to the receiver placement for the estimation of dynamic elastic moduli. Specifically, for Case 1, Receiver 4 happens to be along the line between the excitation and crack center, and directly below the center of the crack as shown in Figure 2, so we could accurately estimate the dynamic elastic moduli. While, for the first arrangement of Case 2, along the path from the excitation to Receiver 4, seismic wave doesn't pass through cracks as shown in Figure 12(a), so the estimated dynamic elastic moduli are approximately consistent with the value of homogeneous media. For the second arrangement of Case 2, the crack centers are arranged along the path from excitation to Receiver 4 in a stacked manner. From the wavefield snapshots as shown in Figure 12(b), stacked cracks prevent Receiver 4 from receiving the Page 11 of 18 SM-DFN modeling for seismic responses (a)L e= (b)L e= �N oqpal�ed�isp 亏 emen� 。 t No� 皂 l�ed�isp 亏 emen� 。。 f Et (c)L e= 即 I Figure 11: Snapshots of the vertical displacement wavefield at t = 0.70 s: LSM-DFN models with crack half-lengths of 50 m, 100 m, 200 m and four different inclined angles of 30 ° , 45 ° , 60 ° , and 90 ° .direct seismic signal excited by the source, so the dynamic elastic modulus cannot be accurately estimated based on this receiver. Figure 13(b) shows the average dynamic elastic moduli obtained from 3 receivers, Receivers 3-5 where the tendency seems more reasonable. Therefore, we may conclude that the obtained dynamic elastic moduli mainly depend on the receiver arrangement and accurate travel-time estimation. Figure 14 and Figure 15 show the differential stress fields under quasi-static loading state at Ey = 0.60%, and snapshots of the vertical displacement wavefield at t = 0. 70 s of the LSM-DFN models with inclined angles of 30 ° , 45 ° ,60 ° , and 90 ° in three types of spatial arrangements. It can be seen from Figure 14 that as the distance between the cracktips decreases, the stress concentration phenomenon increases. From the perspective of the seismic wavefield, Figure 15 further explains the diverse degrees of static and dynamic stress concentration caused by fracture arrangements and inclined angles. Different crack arrangements with various inclined angles make the medium show varied equivalent elastic modulus, resulting in uneven differential stress field distribution, and around the crack tips show distinct static and dynamic stress concentrations. Such numerical simulations provide valuable insights into the mechanisms and processes related with stress variation and seismic responses. Moreover, they could help assess theories by comparing simulations with predictions, serving as powerful tools to improve both theory and experiment. 4.3. Case 3: LSM-DFN modeling for naturally-fractured reservoirs According to Case 1 and Case 2, we examine the validity in emulating elastic wave propagation and scattering in cracked media using LSM-DFN model and analyze the effects of the receiver placement on the dynamic elastic moduli. In this case, we apply the integrated workflow of LSM-DFN modeling proposed by Liu and Fu (2020b) to Page 12 of 18 SM-DFN modeling for seismic responses
N OIJilalieed�isp 亏 ement, 。。 t印 (a)Coplanar cracks N oqpal�ed 4isp�cemen� 。 T 。 .. 厅仆一 (b)Stacked cracks(c)Mixed cracks : :I Figure 12: LSM-DFN model with a set of four cracks with a half-length of 50 m and inclined angle of 0 ° in three types of spatial arrangements, (a) coplanar, (6) stacked, and (c) mixed cracks: Differential stress fields under quasi-static loading state at Ey = 0.60% , and snapshots of the vertical displacement wavefield at t = 0.70 s (from left to right). seismic responses in naturally-fractured reservoirs. Figure 16 shows two real outcrops of Ordovician carbonates in the northwest of the Shuntuoguole low uplift of Tarim Basin in China. The complex fractured reservoirs with different scales and orientations manifest multi-stage tectonic movements (Liu and Fu, 2020b). The complex fractured systems raise challenges in high-accurate modeling for natural fracture networks. Liu and Fu (2020b) extracted manually from a digital photograph by the introduction of Healy, Rizzo, Cornwell, Farrell, Watkins, Timms, Gomez-Rivas and Smith (2017). Then, Liu (2020) improve the fracture extraction by image processing method, the gradient Hough transform (GrdHT) introduced by Mukhopadhyay and Chaudhuri (2015). Figure 17 displays the flow diagrams for efficient line detection, where those fracture networks could be separated from the image background with a gradient magnitude less than a certain threshold. Figure 17(a) shows smaller-scale cracks relative to the overall model size, whilst Figure 17(b) contains multi-scale cracks. In Figure 17(b), those large-scale cracks almost cross the whole sample and divide the sample into nearly intact large blocks. The effective stiffness of the blocky system is seen as stronger than those of samples with homogeneously distributed small fractures (Harthong et al., 2012). Figure 18 shows the distribution characteristics of fracture length of the two outcrops. We get 171 and 619 fractures from the two outcrops, respectively. The fracture lengths of Outcrop 1 are mainly concentrated in a small Page 13 of 18
15 30 45 60 75 900.00.20.40.60.81.01.2
Inclined angle (°)1-s 1-d2-s 2-d3-s 3-dTheoretical solution
Inclined angle (°)1-s 1-d2-s 2-d3-s 3-dTheoretical solution
SM-DFN modeling for seismic responses (a)Coplanar cracks(b)Stacked cracks
Differential I 创 (c)Mixed cracksFigure 14: Differential stress fields under quasi-static loading state at E Y = ° , 45 ° , 60 ° , and 90 ° in three types of spatial arrangements. References
Bonilla-Sierra, V., Scholtes, L., Donze, F., Elmouttie, M., 2015. Rock slope stability analysis using photogrammetric data and dfn-dem modelling. Acta Geotechnica 10, 497-511. Cheng, H., Luding, S., Saitoh, K., Magnanimo, V., 2020. Elastic wave propagation in dry granular media: effects of probing characteristics and stress history. International journal of solids and structures 187, 85-99. Ciccotti, M., Mulargia, F., 2004. Differences between static and dynamic elastic moduli of a typical seismogenic rock. Geophysical Journal International 157, 474--477. Crampin, S., 1978. Seismic-wave propagation through a cracked solid: polarization as a possible dilatancy diagnostic. Geophysical Journal International 53, 467--496. Crampin, S., 1984. Effective anisotropic elastic constants for wave propagation through cracked solids. Geophysical Journal International 76, 135-145.Cundall, P.A., 1971. A computer model for simulating progressive, large-scale movement in blocky rock system, in: Proceedings of the International Symposium on Rock Mechanics, 1971. Cundall, P.A., Strack, O.D., 1979. A discrete numerical model for granular assemblies. geotechnique 29, 47-65. Curtin, W., Scher, H., 1990a. Brittle fracture in disordered materials: A spring network model. Journal of Materials Research 5, 535-553. Curtin, W., Scher, H., 1990b. Mechanics modeling using a spring network. Journal of Materials Research 5, 554-562. Goodman, R.E., Taylor, R.L., Brekke, T.L., 1968. A model for the mechanics of jointed rock. Journal of Soil Mechanics & Foundations Div . Griffiths, D., Mustoe, G.G., 2001. Modelling of elastic continua using a grillage of structural elements based on discrete element concepts. Inter-national Journal for Numerical Methods in Engineering 50, 1759-1775. Harthong, B., Scholtes, L., Donze, F.V., 2012. Strength characterization of rock masses, using a coupled dem-dfn model. Geophysical Journal International 191, 467--480. Hassold, G., Srolovitz, D., 1989. Brittle fracture in materials with random defects. Physical Review B 39, 9273. Healy, D., Rizzo, R.E., Cornwell, D.G., Farrell, N.J., Watkins, H., Timms, N.E., Gomez-Rivas, E., Smith, M., 2017. Fracpaq: A matlab™ toolbox Page 15 of 18
SM-DFN modeling for seismic responses Fl / / // ///1 R / / LJ Fl (a)Coplanar cracks n (b)Stacked cracks I / / / I n f 。 h t g 吆0� 皂 l\8ed�isp 亏 ement-,en 。 f a t 印— h k c a r st c n e a hm e t . wg n a s rr e da 。 a. m t Na p sF fD 。 sM Se LP y .. t se e r 。 h 7 s ck=.In t a ro ct0 da9 e xdd _n a e . f 汕 )e, 。 。 cv a ( w6 t 。 n5 e4 m 。 。 e 3 c a f 。 p s . se d g a cn ·_ t r ed ve ehc tn . f ot n e s r t ofe f h s._pd r a u n 。 s f •• 5d 1n ae urm 。 g ·1 F5 } oqpali.z:ed�isp�cement, 。 TI 讫oqpali.z:ed�isp�cement, 。 T”t田胪田if for the quantification of fracture patterns. Journal of Structural Geology 95, 1-16. Hrennikoff, A., 1941. Solution of problems of elasticity by the framework method. Journal of Applied Mechanics. Huang, Y., Chandra, A., Jiang, Z., Wei, X., Hu, K., 1996. The numerical calculation of two-dimensional effective moduli for microcracked solids. International journal of solids and structures 33, 1575-1586. Hudson, J.A., 1981. Wave speeds and attenuation of elastic waves in material con血ningcracks. Geophysical Journal International 64, 133-150. Hudson, J.A., 1986. A higher order approximation to the wave propagation constants for a cracked solid. Geophysical Journal International 87, 265-274.Hyman, J.D., Gable, C.W., Painter, S.L., Makedonska, N., 2014. Conforming delaunay triangulation of stochastically generated three dimensional 邮cretefracture networks: a feature rejection algorithm for meshing strategy. SIAM Journal on Scientific Computing 36, A1871-A1894. Jiang, C., Zhao, G.F., Khalili, N., 2017. On crack propagation in brittle material using the distinct lattice spring model. International Journal of Solids and Structures 118, 41-57. Jing, L., Hudson, J., 2002. NUillerical methods in rock mechanics. International Journal of Rock Mechanics and Mining Sciences 39, 409--427. Karihaloo, B.L., Shao, P., Xiao, Q., 2003. Lattice modelling of the failure of particle composites. Engineering Fracture Mechanics 70, 2385-2406. Kosteski, L., D'ambra, R.B., lturrioz, I., 2009. Fractomechanics parameter calculus using the discrete element method with bars. Latin American Journal of Solids & Structures 6. Lei, Q., Latham, J.P., Tsang, C.F., 2017. The use of discrete fracture networks for modelling coupled geomechanical and hydrological behaviour of fractured rocks. Computers and Geotechnics 85, 151-176. Lilliu, G., van Mier, J.G., 2003. 3d lattice type fracture model for concrete. Engineering Fracture Mechanics 70, 927-941. Lisjak, A., Grasselli, G., 2014. A review of discrete modeling techniques for fracturing processes in discontinuous rock masses. Journal of Rock Mechanics and Geotechnical Engineering 6, 301-314. Liu, N., 2020. Elastic characteristics for naturally-fractured reservoirs using an integrated lsm-dfn scheme with rough contact deformation. arXiv preprint arXiv:2011.01760. Liu, N., Fu, L.Y., 2020a. Elastic characteristics of digital cores from longmaxi shale using lattice spring models. Communications in Computational Page 16 of 18
SM-DFN modeling for seismic responses (a)Outcrop 1 (b)Outcrop 2Figure 16: Naturally-fractured media of two realistic outcrops from Ordovician carbonates in the northwest of the Shuntuoguole low uplift of Tarim Basin in China. Physics 28, 518-538. Liu, N., Fu, L.Y., 2020b. Stress-orientation effects on the effective elastic anisotropy of complex fractured media using the lsm-dfn model. Interpretation 8, SP31-SP42. Liu, N., Fu, L.Y., Cao, C.H., Liu, J.L., 2020a. Research on numerical modeling method of lsm-rvm and toc content influence for digital core from longmaxi formation shale. Chinese Journal of Geophysics-Chinese Edition 63, 2774-2785. Liu, N., Fu, L.Y., Tang, G., Kong, Y., Xu, X.Y., 2020b. Modified lsm for size-dependent wave propagation: comparison with modified couple stress theory. Acta Mechanica 231, 1285-1304. Monette, L., Anderson, M., 1994. Elastic and fracture properties of the two-dimensional triangular and square lattices. Modelling and Simulation in Materials Science and Engineering 2, 53. Mora, P., Place, D., 1998. Numerical simulation of earthquake faults with gouge: toward a comprehensive explanation for the heat flow paradox. Journal of Geophysical Research: Solid Earth 103, 21067-21089. Mukhopadhyay, P., Chaudhuri, B.B., 2015. A survey of hough transform. Pattern Recognition 48, 993-1010. Murai, Y., 2007. Scattering attenuation, dispersion and reflection of sh waves in two-dimensional elastic media with densely distributed cracks. Geophysical Journal International 168, 211-223. Nayfeh, A.H., Hefzy, M.S., 1978. Continuum modeling of three-dimensional truss-like space structures. AIAA Journal 16, 779-787. O'Brien, G.S., Bean, C.J., 2004. A 3d discrete numerical elastic lattice method for seislnic wave propagation in heterogeneous media with topog-raphy. Geophysical Research Letters 31. Ostoja-Starzewski, M., 2002. Lattice models in micromechanics. Applied Mechanics Reviews 55, 35--60. Otter, J.R.H., Cassell, A.C., Hobbs, R.E., 1966. Dynamic relaxation. Procee山ngsof the Institution of Civil Engineers 35, 633--656. Pan, Z., Ma, R., Wang, D., Chen, A., 2018. A review oflattice type model in fracture mechanics: theory, applications, and perspectives. Engineering Fracture Mechanics 190, 382-409. Pao, Y.H., Mow, C.C., 1971. Diffraction of elastic waves and dynamic stress concentrations. Russak. Saint Venant, B., 1863. Sur la distribution des elasticites autour de chaque point d'un solide ou d'un而lieude contexture quelconque, particul论rement lorsqu'il est amorphe sans etre isotrope. Journal de Math. Pures et Appliquees 8, 257-430. Saltzer, S.D., Pollard, D.D., 1992. Distinct element modeling of structures formed in sedimentary overburden by extensional reactivation of basement normal faults. Tectonics 11, 165-174. Schnaid, F., Spinelli, L.d.F., Iturrioz, I., Rocha, M., 2004. Fracture mechanics in ground improvement design. Proceedings of the Institution of Civil Engineers-Ground Improvement 8, 7-15. Suiker, A., Metrikine, A., De Borst, R., 2001. Comparison of wave propagation characteristics of the cosserat continuum model and corresponding discrete lattice models. International Journal of Solids and Structures 38, 1563-1583. Swoboda, G., Shen, X., Rosas, L., 1998. Damage model for jointed rock mass and its application to tunnelling. Computers and Geotechnics 22, 183—203. Vlastos, S., Liu, E., Main, LG., Narteau, C., 2007. Numerical simulation of wave propagation in 2-d fractured media: scattering attenuation at different stages of the growth of a fracture population. Geophysical Journal International 171, 865-880. Walsh, J., 1965. The effect of cracks on the uniaxial elastic compression of rocks. Journal of Geophysical Research 70, 399-411. Wang, X., Cai, M., 2019. A dfn-dem multi-scale modeling approach for simulating tunnel excavation response in jointed rock masses. Rock Mechanics and Rock Engineering, 1-25. Zhao, G.F., 2010. Development of Micro-Macro Continuum-Discontinuum Coupled Numerical Method. Ph.D. thesis. Ecole Polytechnique FEDErale de Lausanne. Zhao, G.F., 2015. Modelling 3d jointed rock masses using a lattice spring model. International Journal of Rock Mechanics and Mining Sciences 78, 79-90. Page 17 of 18
SM-DFN modeling for seismic responses Trace length map, n =
171 Trace length map, n = 。 \勹 二 I' 、 严 、I � / 一 二
200 1-,..,---/ \.� \� '\:"" --100 。。 S;!ll;!UI'A —---- -v--. 二 勹 < 斗/ _..--,"-义 召 _..) '-\-� 一 - / X-" .. ---....., _ / \ 琴 - 之 / _LJ 三 I ,'e 立-- _,, --、一 又 -------- - <'... ' 一 一 - ------ 、- 一、 /夕 一气 -< /1 - \... 、 _ - - -、文
600 100 玄今 ---- -"\ \ 忑w,-L_ '., 之 ^ \...
500 X, metres 600
500 600 20 40 60 80Trace length, metres (a)Outcrop 1
100 120 50 100 150 200 250 300 350Trace length, metres (b)Outcrop 2
Figure 18: Crack length statistics of (a) Outcrop 1 and (b) Outcrop 2 t: = y i:; Y= ,血.. ,. , .4,`上书\· ``量, .. •• J •• f ;.'{`` , ,`d, 乒 it. `会 }. 、知 ,· ... ` \ 蝎. `”,4' 生,. 、 C ·f,kPri ., -·· c 谧 圈 趴 = (a)Outcrop 1 t = = = > .、 >o'- "'� . - ,, � '½ - -- '• , < � '--o - 夕 、 , I • . ✓ ------ 、 - . , -- , - 只义· -、 - ••• 干 - /< --:-- --', 、 一 , 夕 - 'c. 一了二一 =� .. 矗.' / - r ,- - 雹 � ' 一了 ·- ....... --< -_ ,~ �,'c- ''' 一 < - - / ' / \、,_ -� ' '、 - - 一 千-习�' 圈 -- .. ' ' < � ---- �-- ' .. \ 又--三` - � -· -· , --. - 式 �- -< � -.、-- / � 一 , -仁 ` -、 __. , 俨 - :;,_ 之 _, "-了 —一 气三了二三冷乍 --_ f y= l.00% Ff 『 。乌 "'---0 � I» t = l.20 s 示 Oz i � 旦 了 § 、 了 -7 ""' - ~ - [ ~ 一 -� \ ,,.., 扛令 - 之 /---'t':�_o了三、飞尸 (b)Outcrop 2 Figure 19: Numerical results of Outcrop 1: (a) Differential stress fields under quasi-static loading state at Ey = 0.25%, 鸟= 0.50%, Ey = 0.75%, and Ey = 1.00%; (b) Snapshots of the vertical displacement wavefield at t = 0.60 s, t = 0.80 s, t = 1.00 s, and t = 1.20 s. Page 19 of 18 .0 0.2 0.4 0.6 0.8 1.0-50510152025 Normalized timeSource Receiver 1 Receiver 2Receiver 3 Receiver 4 Receiver 5Receiver 6 Receiver 7 Receiver 80.0 0.2 0.4 0.6 0.8 1.0-50510152025