Displacement Discontinuity Research Articles

A bulking model to simulate stope convergence in deep tabular excavations

This paper investigates the use of a novel time-dependent bulking model in a displacement discontinuity code to simulate the convergence between the roof and the floor of deep tabular stopes. This study is important as the current layout design criteria that are based on elastic theory have become outdated. This approach ignores rock failure and fails to predict the actual magnitude and the time-dependent nature of the stope convergence. This has implications for support design and for estimating the stress acting on remnants in mature mines. A rock mass modulus reduction strategy, often used to circumvent this problem, reflects the expected lower effective rock modulus for near-stope fractured material, but must be applied for the entire rock mass region in the displacement discontinuity solution procedure for tabular layout problems. As an alternative, it is proposed in this paper that the rock bulking can be expressed as a function of the elastic closure by using a constitutive rule for the reef-normal compressive reaction stress that is expressed as a function of the compaction strain. The proposed model was tested using convergence data collected from a shaft pillar extraction in a deep gold mine. The bulking model, in conjunction with a limit equilibrium model to simulate the face crushing, seems capable of replicating the underground behaviour. As the stresses on remnants and pillars can be more accurately simulated, this approach holds promise of improving current design criteria. Calibration of the model nevertheless remains a challenge and a programme of routine convergence measurements will have to be implemented on the mines.

A computational framework for crack propagation along contact interfaces and surfaces under load

We present the first implicit computational framework for simulating crack propagation along contact interfaces and surfaces under load in three-dimensional bodies, which is distinct from modelling the contact interaction associated with crack closure. We restrict ourselves to brittle fracture and frictionless contact and focus on numerical challenges associated with the coupling of unilateral constraints emerging from the Griffith’s criterion and the contact conditions. The formulation is based on the configurational mechanics framework and is solved using the finite element method. The approach utilises a monolithic Arbitrary Lagrangian–Eulerian formulation permitting simultaneous resolution of crack propagation and unilateral contact constraints. Contact is embedded in the model using the well-known mortar contact formulation. Evolving cracks are explicitly modelled as displacement discontinuities within the mesh. Heterogeneous approximation of arbitrary order is used to discretise spatial displacements, enabling hp-adaptive refinement around the crack front and the contact interfaces traversed by the crack. The result is a holistic approach which handles issues associated with thermodynamic consistency, numerical accuracy and robustness of the computational scheme. Several numerical examples are presented to verify the model formulation and implementation; they also highlight how contact pressure and load applied on surfaces traversed by cracks influence their propagation. The robustness of the approach is validated by comparison of our simulations with existing numerical results and an industrial experiment involving cracks of complex morphologies propagating along contact interfaces between multiple deformable bodies.

Dynamic stability assessment of consequent rock slopes considering multiple transmissions and reflections of stress waves

In this work, Based on the conservation of momentum at the wave fronts and the displacement discontinuity method, a theoretical approach for calculating the dynamic safety factor of consequent rock slopes is proposed. The multiple transmissions and reflections of stress waves on the joints and free surfaces are taken into account. A comparison of UDEC simulation and theoretical calculation was carried out to verify the feasibility of the proposed method. Moreover, parametric studies were conducted to investigate the effects of free surfaces, joint stiffness, amplitude of incident waves, and frequency of incident waves. The results show that the proposed method can predict the time history curve of the safety factor of a slope readily and accurately. Waves reflected from both the slope and top surfaces are found to have a significant effect on the stability of the slope – considering them or not can cause a difference in the safety factors up to 32.5%. The greater the joint stiffness is, the greater the safety factor fluctuates with time. However, when the joint stiffness is sufficiently high, the effect of joint stiffness on the safety factor can be neglected. The approach proposed in this work provides a theoretical tool for the dynamic assessment and design of consequent rock slopes. Highlights An analytical method is proposed for calculating the dynamic safety factor of rock slopes. Multiple transmissions and reflections of stress waves are considered. The proposed method has high calculation accuracy and can be applied in practice.

Numerical modeling of low-frequency distributed acoustic sensing signals for mixed-mode reactivation

In recent years, the development of distributed acoustic sensing (DAS) technology has enabled direct monitoring of subsurface strain during hydraulic fracturing (HF) operations. Most low-frequency (<1 Hz) DAS (LFDAS) signals exhibit strain or strain-rate patterns that are characteristic of propagating tensile hydraulic fractures. However, mixed-mode fault reactivation, consisting of shear slip (mode II) with tensile opening (mode I), can occur in cases where propagating hydraulic fractures intersect preexisting natural fractures or faults. In this study, we show an anomalous LFDAS signal that was observed during an HF operation, with characteristics that differ from typical signals from tensile hydraulic fractures. Anomalous characteristics include the onset of an asymmetric compression-extension doublet after pumping was terminated. The location of the signals, coupled with the image log and seismic data, suggests that mixed-mode reactivation occurred on a preexisting fault. We use a simplified numerical model based on the displacement discontinuity method (DDM) to simulate the anomalous LFDAS signal. Results find that the first-order characteristics of the anomalous signal can be approximated as an initial tensile fault opening (mode I) followed by a shear slip (mode II) on a fault. Therefore, we demostrate the approach of using DDM to investigate mixed-mode fault reactivation during HF operations.

A concise review of small-strain phase-field modeling of ductile fracture

The computational modeling of ductile fracture is a complex problem due to the simultaneous presence of various competing dissipation mechanisms at the microscale resulting in evolving displacement discontinuities at the macroscale. Phase-field modeling of brittle fracture has emerged as one of the more computationally efficient smeared approaches for simulating crack propagation in solids. In recent years, there has been significant attention towards extending this approach to ductile fracture. The first papers on this topic were published approximately a decade ago, and since then, there has been a growing number of new contributions proposing different approaches and applications. The purpose of this work is to offer an extensive bibliography on phase-field modeling of ductile fracture, together with a concise presentation of some of the most widely used formulations, specifically limited to small-strain plasticity under monotonic loading. Few aspects of particular interest are also discussed in detail.

Numerical simulation of proppant migration in horizontal wells with multi-fracture fracturing

For horizontal wells fracturing in unconventional reservoirs, the fracture morphology and proppant migration are severely affected by the interfering effects among multiple fractures propagation. Current numerical simulation for proppant migration did not reveal the interrelation among flow distribution, fracture length and sandbank shape during multiple fractures propagation. In this study, a proppant migration model is proposed by combining the displacement discontinuity method (DDM) and the Eulerian-Lagrangian method (E-L M), which implement dynamic multiple fractures propagation simulation as well as the changing state of flow distribution for each fracture. The accuracy of the numerical model is verified by laboratory experiments. The sensitivity parameters of proppant settlement are analyzed, and two dimensionless parameters, height ratio γ and area ratio η are proposed for proppant migration evaluation. The results show that with a lower fracturing fluid flow rate, fluid viscosity, the height ratio γ is lower, and the area ratio η is higher, which also results in the proppant being more likely to accumulate at the fracture entrance. As the proppant particle size increases from 100 μm to 500 μm, the dimensionless height ratio γ and the dimensionless area ratio η decrease significantly. Increase in fracture spacing during horizontal well multi-cluster fracturing results in a limited improvement for proppant settlement. Formation elastic modulus has an appreciable impact on proppant migration in multi-fracture propagation, high elastic modulus formation leads to an ununiform proppant placement over long distances.

Displacement Discontinuity Method Taking into Account the Curvature of the Crack

Displacement Discontinuity Method Taking into Account the Curvature of the Crack

An unequal fracturing stage spacing optimization model for hydraulic fracturing that considers cementing interface integrity

Determining reasonable fracturing stage spacing is the key to horizontal well fracturing. Different from traditional stage spacing optimization methods based on the principle of maximum stimulated reservoir volume, in this paper, by considering the integrity of the wellbore interface, a fracture propagation model was established based on displacement discontinuity method and the competition mechanism of multifracture joint expansion, leading to the proposal of an unequal stage spacing optimization model. The results show that in the first stage, the interfacial fractures spread symmetrically along the axis of the central point during that stage, while in the second and subsequent stages, the interfacial fractures of each cluster extend asymmetrically along the left and right sides. There are two kinds of interface connectivity behaviour: in one, the existing fractures first extend and connect within the stage, and in the other, the fractures first extend in the direction close to the previous stage, with the specific behaviour depending on the combined effect of stress shadow and flow competition during hydraulic fracture expansion. The stage spacing is positively correlated with the number of fractures and Young’s modulus of the cement and formation and is negatively correlated with the cluster spacing and horizontal principal stress difference. The sensitivity is the strongest when the Young’s modulus of the cement sheath is 10–20 GPa, and the sensitivity of the horizontal principal stress difference is the weakest.

Study of stress wave propagation across a non-persistent joint based on a boundary integral equation method

Studies on stress wave propagation across persistent joints have been conducted extensively. Nevertheless, there exists a consensus that non-persistent joints are widely and densely distributed, which have a profound impact on wave propagation in jointed rock masses. A boundary integral equation method is suggested in this paper to investigate the characteristics of transmitted wave field for the case of stress wave propagation across a single non-persistent joint. The displacement continuity and discontinuity boundaries are combined in the method. The method presented in the current study is applicable to the analysis of wave propagation across non-persistent joints with arbitrary incident angles. Then, taking a single non-persistent joint arranged with only one joint segment as an example, the applicability of the method in dealing with the problems of wave propagation is verified by comparing the results with those from the discrete element method and analytical methods. Subsequently, parametric studies are carried out, including the effects of joint-segment length, rock-bridge length, wave frequency and incident angle on the transmitted wave. The result indicates that the existence of non-persistent joint makes the transmitted displacement field different from that of persistent joint, because the scattered wave is produced during the process of wave propagation. The displacement amplitude may be amplified evidently in some regions and the spatial distribution pattern of the transmission coefficient is closely related to the joint-segment length, rock-bridge length and incident wavelength.

New crack front enrichment for XFEM modeling

The extended finite element method (XFEM) using a new crack front enrichment is proposed for the crack modeling in this study. The discontinuity of displacement in the element including the crack tip or front is modeled by combining the Heaviside function with a step function. The proposed enrichment scheme avoids the calculation of the geometrical characteristic lengths in the element containing the crack tip or front when using the Ramp function, and the step function only depending on the level set function is defined to activate the discontinuity. Thus, both two dimensional (2D) and three dimensional (3D) cracks can be modeled more effectively by the new crack front enrichment. The Conjunction Gradient Method (CGM) that can reduce the burden of storage and manipulation is utilized to solve the global linear algebraic system of 3D XFEM. In addition, two kinds of 3D variable-node elements are used to connect the coarse and fine elements to save the number of nodes and elements. The Westergaard’s crack problem is analyzed to test the convergence rate of the proposed method, and it is found that the new crack front enrichment can match the theoretical convergence rate of the finite element method (FEM) in the presence of a dominant r-singularity. The proposed method is applied to determine the stress intensity factors (SIFs) of straight, inclined, edge and curved cracks. Comparing with the reference solutions demonstrates that the proposed crack front enrichment exhibits significant advantages for the 3D crack modeling.

Characterization of landslide displacements in an active fault zone in Northwest China

AbstractLandslides can be caused by natural forcing and anthropogenic activities. Zhouqu County (China) on the eastern margin of Qinghai‐Tibet Plateau is set within the active Pingding‐Huama fault zone with evident fractures on the land surface. Frequent landslides and debris flows have occurred in this region due to river erosion, rainfall and deforestation. Here we quantified the slope movements using time‐series synthetic aperture radar interferometry (InSAR) based on the ascending and descending Sentinel‐1 satellite images acquired between October 2014 and August 2020. We observed distinct displacements in the highly fractured fault zone. The eastward and vertical displacement time series between February 2017 and July 2020 were constrained by the common‐day ascending and descending acquisitions. The eastward rates (461 mm/year) were greater than those in the vertical direction (−185 mm/year). We also note displacement discontinuities across the thrust faults beneath the Suoertou and Zhongpai landslides. Seasonal variations in the displacement time series suggest that the cyclic rainfall is the primary driver for the mass wasting processes rather than the tectonic loading. As a complement to in situ observations, our results demonstrate that InSAR is an effective tool to characterize the spatio‐temporal nature of landslide displacements in complicated geological environments.Plain Language SummaryZhouqu County in the Pingding‐Huama fault zone in the eastern margin of Qinghai‐Tibet Plateau is identified as a high priority site to research on clusters of landslides and debris flows in a mixed geodynamic setting of active tectonics, seasonal rainfall, river erosion and anthropogenic activities. However, our knowledge about landslide kinematics in this complicated region is still limited. We relied on remote sensing images from one ascending and one descending Sentinel‐1 satellite tracks to constrain the spatial–temporal displacement dynamics of active landslides from 2014 to 2020. The spatial patterns of displacements are determined by thrust faulting, river erosion, and anthropogenic activities. The temporal variations of landslide speed are mainly controlled by the seasonal rainfall rather than the tectonic loading.

Low-frequency distributed acoustic sensing shape factors for fracture front detection

Accurate knowledge of fracture extents generated in multistage unconventional completions remains elusive. Crosswell low-frequency distributed acoustic sensing (LF-DAS) measurements can determine the time and location of a frac hit. Knowing where and when a frac hit occurs constrains the fracture extent but does not estimate it quantitatively. A recent study on crosswell LF-DAS demonstrated a simple method to rapidly determine the instantaneous fracture propagation rate when a frac hit occurs. This method, the zero strain-rate location method (ZSRLM), is based on laboratory experiments and numerical modeling assuming a radial fracture geometry. The method estimates a fracture propagation velocity that is used to extrapolate the final fracture extent. The propagation rate is calculated based on dynamic estimates of the nearest distance from the fiber to the front of a propagating fracture. The ZSRLM is adapted to estimate the distance to the fracture front based on rectangular fracture geometries. A 3D displacement discontinuity method program generates crosswell LF-DAS strain-rate waterfall plots considering a single, rectangular fracture of constant height. Over 150 different simulations were conducted varying formation mechanical properties, fracture height, and the vertical and horizontal offset between the treatment and monitor well. For each simulated case, the ZSRLM is used to estimate the distance to the fracture front based on the simulated waterfall plots. The difference between the estimated and actual distance to the front is corrected by a shape factor. The relationship among the shape factor, fracture height ratio, and vertical offset ratio is determined. Using a shape factor improves the performance of the ZSRLM by up to a factor of two for rectangular fractures. The updated ZSRLM is applied to extrapolate final fracture extents in two field cases: a single cluster stage in the Montney Formation and a multicluster stage of an Austin Chalk completion.

Analysis of 3D quasi-brittle solids failures by crack growth using the strong discontinuity approach with the boundary element method

The continuum strong discontinuity approach (CSDA) is suitable for the analysis of material failures of quasi-brittle solids when constitutive damage models are applied to kinematic equations with displacement discontinuities. Its association with the implicit version of the boundary element method (BEM), using cells with embedded strong discontinuity, allows the analysis of material failures with crack growth. Such cells are mandatory for the domain discretization, which are restricted to the region where the crack is supposed to occur. They are introduced in the model as the analysis progresses, without the need to re-evaluate the pre-existing terms of the matrices, but only through their appropriate expansion. Discontinuity in the displacement field is quantified, after internal equilibrium verification, assuming a uniform displacement jump inside each cell. This work uses this methodology for the analysis of three-dimensional solids under external loads. An elastic-degrading isotropic constitutive model with an exponential softening law is considered, with the strong discontinuity regime imposed directly after the end of the elastic regime. This is the first time that such cells are applied for full three-dimensional mixed mode fracture problems, in the context of the boundary element method. The results presented in the numerical examples are in good agreement with the references.

Three-dimensional interfacial mechanical-thermal fracture analysis of a two-dimensional hexagonal quasicrystal coating structure

Based on the displacement and temperature discontinuity boundary element method, three-dimensional thermal fracture behaviors of two-dimensional hexagonal quasicrystal (QC) coating with an arbitrary interface crack are studied. First, fundamental solutions of discontinuities are derived using the general solution and Hankel transform technique, and boundary integral-differential equations are constructed at interface. Then, taken a penny-shaped interface crack dispersed by ring elements as an example, Green’s functions with discontinuities are derived on these elements. After that, stress intensity factors (SIFs) and the energy release rate (ERR) are expressed with discontinuities. Finally, effects on fracture were discussed, including coating thickness, applied loads, and material mismatch.

Anelastic dispersion and attenuation of P- and SV-wave scattering by aligned nonisothermal inclusions of finite thickness

We have investigated the anelastic dispersion and attenuation of P- and SV-wave scattering by nonisothermal inclusions of finite thickness. The inclusions, which are aligned and sparsely embedded in an isotropic medium, induce an initial static stress field (acoustoelasticity) and a nonlinear dependence of the velocities on this stress. Moreover, we describe the anelastic properties as a function of frequency by incorporating the displacement discontinuities across the inclusion into the representation theorem by using the Foldy approximation. The response as a function of temperature is calculated for different incidence angles (anisotropy), and the results find that anelasticity increases with an increasing temperature difference between the inclusions and the background medium. The SV wave in the solid inclusions indicates a stronger sensitivity along the inclusion normal, and it is more affected than the P wave. The P- and SV-wave scattering by fluid-saturated inclusions behave in an opposite manner. This theory can be useful to evaluate the distribution of temperature from seismic attributes.

Integrated simulation of three-dimensional hydraulic fracture propagation and Lagrangian proppant transport in multilayered reservoirs

The numerical simulation of hydraulic fractures is a critical yet challenging computational problem due to its multi-physics nature. This paper develops an integrated hydraulic fracturing simulator by coupling a planar three-dimensional (PL3D) fracture propagation model with an efficient Eulerian–Lagrangian (E–L) proppant transport model. The fracture propagation model uses the finite volume method (FVM) and displacement discontinuity method (DDM) to solve the fluid flow and rock deformation, respectively. For proppant transport, we develop a pseudo-3D multiphase particle-in-cell (P3D MP-PIC) model, in which the fluid flow is addressed by FVM as a continuum, but the particles are tracked in a Lagrangian fashion as discrete phases. In contrast to Eulerian–Eulerian (E–E) proppant transport models, the P3D MP-PIC can efficiently deal with multi-modal particle simulations (i.e., particles of different sizes or materials) and avoid the problem of numerical diffusion. The fracture propagation and proppant transport models are validated by analytical solutions and laboratory experiments. We adopt a one-way coupling strategy to consider the effect of complex fracture propagation and fluid leak-off on slurry transport, in which the dynamic fracture geometry and fluid leak-off is first computed via the fracture propagation model, followed by the fully coupled fluid–particle simulation using the P3D MP-PIC. The integrated model can simulate the fracturing treatments in multilayered reservoirs with varying confining stresses at an industrial field scale. The simulation results can improve the prediction of effective/propped fracture geometries and better the fracturing designs.

Geometrically exact beam theory with embedded strong discontinuities for the modeling of failure in structures. Part I: Formulation and finite element implementation

We introduce a three-dimensional geometrically nonlinear Reissner beam theory with embedded strong discontinuities for the modeling of failure in structures and discuss its finite element implementation. Existing embedded beam theories are geometrically linear or two-dimensional, motivating the need for the present work. We propose a geometrically nonlinear beam theory that accounts for cracks through displacement discontinuities in 3D. To represent the three modes of fracture inside an element, we enrich each displacement field with an incompatible mode parameter in the form of a jump, an additional degree of freedom. We then eliminate these nodeless degrees of freedom through static condensation and evaluate them within the framework of inelasticity by utilizing an operator split solution procedure. Seeing that the coupled strain-softening problem is non-convex, we present an alternating minimization (staggered) algorithm, thus retaining a positive definite stiffness matrix. Finally, the four-parameter representation by quaternions describes a three-dimensional finite rotation. We demonstrate a very satisfying and robust performance of these new finite elements in several numerical examples, including the fracture of random lattice structures with application to fibrous materials. We show that accounting for geometrical nonlinearity in the beam formulation is necessary for direct numerical simulations of fiber networks regardless of the density.

Evaluation of Parent Well Depletion Effects on Fracture Geometry Based on Low-Frequency Distributed Acoustic Sensing in Hydraulic Fracture Test Site-2

Summary Low-frequency distributed acoustic sensing (LF-DAS) has recently received much attention for its ability to monitor fracture propagation at offset wells. Hydraulic Fracture Test Site-2 (HFTS-2), a Department of Energy–sponsored field-based research experiment, has acquired LF-DAS data sets during the stimulation of many horizontal wells in the Wolfcamp Formation of the Permian Basin. Over 100 stimulated stages with different completion designs in four horizontal wells were monitored by two horizontal offset wells and one vertical pilot hole with permanent fibers. The parent well depletion affected all four horizontal wells in almost half of their lateral section. Several studies have been performed on the acquired comprehensive data set. In this study, we apply our novel Green’s function-based inversion algorithm to calculate the fracture width of each stage and investigate the impact of parent well depletion on fracture geometry. The Green’s function-based inversion algorithm has been successfully applied to a few stages of field case studies. The Green’s function matrix was built based on linear elasticity theory (3D displacement discontinuity method) to relate the fracture openings with strain responses measured along the length of the fiber. This novel algorithm only relies on measured strain to solve fracture geometry at the monitoring well. Therefore, it is independent of the physics of the fracture propagation process and can be used to validate hydraulic fracture modeling results. Using our inversion algorithm, we can efficiently and quantitatively interpret LF-DAS data to provide information on fracture geometry and completion efficiency. We analyze more than 100 stages of the LF-DAS measurements obtained at the fiber wells B3H and B4H during B1H, B2H, and B4H stimulations. We apply our inversion algorithm for four data sets covering the stimulation of the abovementioned three wells. The fracture growth at stages in the parent well’s depleted zone is biased more toward upper formation and in the easterly direction. The fracture width at the stages in the parent well’s depletion zone is reduced compared to fracture widths at stages in the nondepleted zone irrespective of the monitoring well location relative to the treatment wells. This difference in fracture widths will affect the proppant distribution to a great extent, thereby affecting the effectiveness of the stimulation. We also illustrate the application of the inversion algorithm for stages that have both conventional fracture hits with “heart-shaped” signals as well as fracture reopening signals. Our inversion algorithm gives a reasonable estimate of the fracture width that aligns with the qualitative analysis of microseismic data sets and statistic summary of fracture hit numbers of HFTS-2. We believe that this quantitative study provides us with insights into the fracture geometry due to parent well depletion effects.

Numerical simulation of fatigue crack propagation in heterogeneous geomaterials under varied loads using displacement discontinuity method

Numerical simulation of fatigue crack propagation in heterogeneous geomaterials under varied loads using displacement discontinuity method

Real-Time Simulation of Hydraulic Fracturing Using a Combined Integrated Finite Difference and Discontinuous Displacement Method: Numerical Algorithm and Field Applications

Real-time simulation of hydraulic fracturing operations is of critical importance to the field-scale stimulation applications. In this paper, we present an efficient yet reasonably accurate program for the numerical modeling of dynamic fractures. Our program, named as FracCSM, is based on combined Integrated Finite Difference (IFD) method and Discontinuous Displacement Method (DDM). FracCSM simulates the initiation and propagation of hydraulic fractures with DDM and mass/heat transport inside fractures by IFD. The frictional loss within the wellbore is also taken into consideration. In this way, we are able to model the propped height and length of the fractures subject to the stress interference effect. Moreover, FracCSM captures the stress shadow effect of multi-stage fractures. To facilitate the monitoring and decision making during the hydraulic fracturing process, we have developed a general framework that supports real-time simulation of fracture propagation. Our developed program demonstrates sound accuracy in comparison with existing simulators. The novelty of this work is the combined simulation algorithm to simulate the multiphysical process during hydraulic fracturing operations. We will demonstrate the program structure as well as the field applications of FracCSM to the real-time simulation of hydraulic fracturing operations in Sulige tight sandstone reservoir.